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NUMERICAL STUDY OF THREE-DIMENSIONAL SEPARATION 
AND FLOW CONTROL AT A WING/BODY JUNCTION 


ABSTRACT 

The problem of three-dimensional separation and flow control at a 
wing/body junction has been investigated numerically using a three- 
dimensional Navier-Stokes code. The numerical code employs an algebra- 
ic grid generation technique for generating the grid for unmodified 
junction and an elliptic grid generation technique for filleted fin 
junction. The results for laminar flow past a blunt fin/flat plate 
junction demonstrate that after grid refinement, the computations agree 
with experiment and reveal a strong dependency of the number of vorti- 
ces at the junction on Mach number and Reynolds number. The numerical 
results for pressure distribution, particle paths and limiting stream- 
lines for turbulent flow past a swept fin show a decrease in the peak 
pressure and in the extent of the separated flow region compared to the 
laminar case. The results for a filleted juncture indicate that the 
streamline patterns lose much of their vortical character with proper 
filleting. Fillets with a radius of three and one-half times the fin 
leading edge diameter or two times the incoming boundary layer thick- 
ness, significantly weaken the usual "necklace" interaction vortex for 
the Mach number and Reynolds number considered in the present study. 
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Chapter 1 
INTRODUCTION 

Three-dimensional separation at intersecting surfaces such as 
wing-fuselage and wing-pylon junctions, and air breathing engine 
inlets is a complex physical phenomenon. The intersecting surface 
flow exhibits complex flow patterns which may significantly affect the 
performance of the aircraft. Major effects include increased pressure 
and heat transfer in the vicinity of the intersections and loss of 
control effectiveness due to flow separation. A variety of simplified 
geometrical configurations have been studied to gain insight into the 
physical aspects of three-dimensional separated flows including (1) 
sharp fin mounted on a flat plate, (2) blunt fin mounted on a flat 
plate, (3) swept compression corner, and (4) semicone affixed to a 
flat plate. 

In the present work we concentrate our attention to the flow past 
a blunt fin/flat plate junction whch is a typical example of junction 
flow. When a high-speed flow passes over a blunt obstruction mounted 
on a flat plate the detached bow shock formed ahead induces a strong 
adverse pressure gradient to the plate boundary layer. This 
disturbance propagates through the subsonic portion of the boundary 
layer and can result in a 3-D separated flow region composed of 
horseshoe vortices near the plate surface and a lambda-type shock 
pattern ahead of the fin. These vortices entrain a part of the free- 
stream flow and spiral downstream. The separation shock emanating 
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from the separated flow region impinges on the fin shock producing a 
shock-shock interaction. Depending on the severity of the incoming 
shock, the shock-shock interaction leads to the formation of 
supersonic jet. This interaction leads to intense heating and high 
pressure locally around the fin leading edge. 

Extensive experimental and computational studies pertinent to the 
physics of the problem have been carried out in the past for flow past 
a wing/body junction and these are reviewed briefly in Chap. 2. The 
description of the governing equations and all possible boundary 
conditions employed in the study is given in Chap. 3. The development 
of the finite volume concept, method of solution and operation of the 
numerical code are given in detail in Chap. 4. The specific physical 
conditions of investigation are described in Chap. 5 and discussion of 
all results are provided in Chap. 6. 



Chapter 2 

PHYSICS OF WING/BODY JUNCTURE FLOWS 
In this chapter, the physics of wing/body juncture flows at low 
speed and high speed are discussed. Key findings from the previous 
studies have been included in the physics of the problem to indicate 
the characteristics of the junction flow. The discussion begins with 
detailed description of low speed unmodified juncture flows for lami- 
nar, turbulent and transitional conditions. Next, the description of 
controlled juncture flows are presented. The high speed unmodified 
and modified juncture flows is considered next to clearly highlight 
the difference in the flow structure with reference to low speed. 
Finally the status of turbulence modeling as applied to the compu- 
tation of junction flow is briefly reviewed. 

2.1 Low Speed Flows 

2.1.1 Low Speed Unmodified Juncture Flows 

If a protuberance is placed on a flat plate over which there is a 
fluid flow, the protuberance produces an adverse pressure gradient in 
its vicinity. This adverse pressure gradient created upstream of the 
protuberance causes the plate boundary layer to deviate from the 
surface resulting in a three-dimensional separation. The separated 
boundary layer rolls up downstream of the separation line to form a 
system of vortices. These vortices are swept around the base of the 
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protuberance and assumes a characteristic shape which has led to its 
name the horseshoe vortex. 

Norman [1]* investigated experimentally the basic features of 
secondary flows and the conditions under which these flows led to 
transition and turbulence. The experimental study defined and des- 
cribed qualitatively the vortex system, the horseshoe vortex system, 
the spiral ing-trai 1 ing vortex system, the vortex system resulting from 
the instability of the top shear layer, and the interactions of these 
sytems on transition. Schwind [2] studied the low speed junction flow 
using smoke visualization technique and observed five types of vortex 
regimes. Schwind observed that as the velocity is increased the 
number of vortices observed also increases and at higher velocities 
the two clockwise rotating vortices began to oscillate with amplitude 
increasing with velocity. Schwind further points out that although 
the flow passed from regime 1 to regime 5 types as the velocity is 
increased, the velocity at which the flow changed from one regime to 
another was by no means well defined and different flow regimes could 
exist at one velocity. Baker [3] demonstrated experimentally the 
existence of both steady and unsteady vortex systems. The oscillatory 
behavior of the horseshoe vortex was found to be complex. As the 
Reynolds minber is increased the steady horseshoe vortex system Degan 
to oscillate intermittently and randomly at two different frequencies. 
As the Reynolds number is increased further, the period of oscillation 
of the vortices became longer, the high frequency oscillation becoming 
more common. At yet higher Reynolds number the oscillations became 
irregular and the horseshoe vortex system became turbulent. It was 
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concluded that the frequency of these oscillations is solely deter- 
mined by the value of Reynolds number and the ratio of diameter of the 
protuberance to the incoming boundary layer thickness. 

Kaul et al. [4] studied numerically the incompressible junction 
flow using the Navier-Stokes equations. The equations were cast in 
generalized coordinates, and solved in time as a hyperbolic system by 
adding a pressure term in the continuity equation. The equations were 
marched until steady state was reached. The calculations indicated a 
strong non-uniform pressure loading along the length of the cylinder. 

A new mechanism for the existence of recirculation bubbles behind the 
cylinder end walls with relatively low ratio of cylinder height to the 
incoming boundary layer thickness was observed. 

The turbulent flow around a wing/body junction is much more 
complex than the laminar counterpart. It consists of a turbulent 
boundary layer interacting with a viscous wake that is pressure driven 
in some regions, Reynold's stress driven in other regions and viscous 
dominated in most other regions. The essential features of the flow 
and the horseshoe root vortex are shown in Fig. 2.1. The wing is 
partially submerged in the oncoming flat plate boundary layer. The 
body is sufficiently long that a turbulent boundary layer exists up- 
stream of the wing. The flow around the wing outside the boundary 
layer is assumed to be inviscid. The vorticity present in the incom- 
ing boundary layer wraps around the wing leading edge to form a horse- 
shoe vortex. The oncoming fluid is decelerated by the adverse pres- 
sure gradient created due to the presence of the wing leading to 
three-dimensional separation. It is in this part that the primary 
root vortex, and possibly secondary vortices are established. 
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Fig. 2.1 Flow structure at the low speed wing/body junction. 
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Meanwhile, fluid off to the side of the wing, passes without 
experiencing much disturbance. Consequently, the fluid elements in 
one part of a vortical line are decelerated relative to fluid in other 
parts. The subsequent wrapping of the vortical line elements around 
the wing leading edge results in the appearance of horseshoe vortex. 
The effect of this skewing of the transverse vorticity line results in 
the generation of streamwise vorticity and a three-dimenstional flow 
velocity. Whereas, the flow upstream of the wing has one component 
velocity, which is directed in the streamwise direction, the flow in 
the root vortex has components in the plane normal to the streamwise 
direction due to the presence of streamwise vorticity. The horseshoe 
vortex experiences viscous diffusion and the effects of counter 
rotation takes place because of the anisotrophy of the Reynold's 
stresses. The flow slowly varies in the streamwise direction and is 
in some sense in equilibrium. The flow does not return to its 
original state downstream from the wing. Further, there are complex 
interactions between the vortices, the separated flow near the trail- 
ing edge and the viscous wake of the wing. The vortices remain 
immersed in the boundary layer from which they originated, and they do 

I 

not lift up away from the body. 

Shabaka and Bradshaw [5, 6] and Oguz [7] made mean flow and 
turbulence measurements in wing-fuselage junctions. The special 
emphasis of the experiment was to investigate the relative importance 
of turbulence terms and examine the spatial extent of the interaction 
flow field. McMahon et al. [8] studied the details of turbulent flow 
in a wing-fuselage type junction using surface oil-flow visualization 
technique combined with hot wire anemonemtry. The experimental data 
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collected includes measurement of three mean velocity components and 
six Reynolds stresses at selected streamwise locations. To complete 
the entire picture of the flow field some measurements were taken 
inside the reversed flow region using hot-wire techniques. It was 
shown that unambiguous measurements of the mean flow direction can be 
made in three-dimensional separated flows. The collected data have 
been used to deduce the flow behavior in and around the juncture. 

Also, the similarities between the mean-flow strain rates and the 
turbulent stresses, and the relative trends among the turbulence 
quantities have been discussed. Hsing and Teng [9] carried out 
experiments on the behavior of three-dimensional turbulent boundary 
layer in the simplified wing-body junction at low subsonic speed. The 
results of the boundary layer behavior in the attached flow region 
have revealed the effects of pressure gradient and curvature of the 
streamline. Briley and McDonald [10] studied numerically the wing/ 
body junction flow with swept leading edges using the Navier-Stokes 
equations. A general nonorthogonal coordinate system that resolved 
the viscous sublayer was employed. The nunerical predictions did not 
agree with measurements well downstream of the leading edge. The 
algebraic turbulence model employed in the study was fairly crude, but 
false diffusion errors, notoriously large in the recirculating region 
may have also contributed. 

2.1.2 Low Speed Controlled Juncture Flows 

The physics of low speed controlled juncture flows differ slight- 
ly from unmodified juncture flows in that the separation at the wing 
leading edge is eliminated by suitable geometrical modifications to 
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the wing. Such modifications include fillets, sweep and suction holes 
in the wing. These simple techniques help in improving the air flow 
around the wing/fuselage junctions and can be easily incorporated in 
existing fleets of aircraft. Preliminary results [11] indicate that 
there is some reduction in juncture drag at moderate angles of attack 
with the use of leading edge fillets. The use of sweep and suction at 
the wing leading edge considerably reduces the extent of separation 
and helps to maintain laminar flow around the juncture. 

Mehta [12] experimentally investigated the effect of wing nose 
shapes on the nature of the generated secondary flow. The 
experimental study investigated in detail the correlation of the 
strength, size and position of the horseshoe vortex with the nose 
shape of the wing. The study clearly demonstrated that the nature of 
the horseshoe vortex system in a wing/body junction is strongly 
dependent on the wing nose shape with size and strength of the vortex 
system increasing with the nose bluntness. In most practical cases, 
this vortex flow persists right up to the wing trailing edge, with 
very little attenuation, so the only practical way to control it is by 
suitably adjusting the wing leading edge shape. Kubendran and Harvey 
[11] have made measurements in the wake region of a simulated wing- 
fuselage juncture, with and without leading-edge fillets, to 
understand the effect of leading edge modifications on the flow field 
around the juncture. Preliminary results have indicated some 
reduction in juncture drag using leading edge fillets at moderate 
angles of attack. Some improvements in the flow characteristics on 
the surface downstream of the juncture were also noticed when leading 
edge fillets were used. As the fillet size was increased, flow 
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characteristics start deteriorating at some point, and an optimum 
fillet size was needed to achieve an overall improvement in the flow 
field. 

The numerical study of controlled juncture flows using fairings 
at the junction have been studied by Sung and Lin [13] and Kubendran 
et al. [14] using the incompressible Navier-Stokes equations. It has 
been found that the leading edge fairing is quite effective in reduc- 
ing the leading edge flow separation and non-uniformity of the wave 
velocity profile and the addition of the trailing edge fairing can 
contribute to further improvement. 

2.2 High Speed Flows 

2.2.1 High Speed Unmodified Juncture Flows 

High speed flow past a blunt fin mounted on a flat plate, results 
in a complex, three-dimensional, viscous-inviscid interaction flow 
field. The detached bow shock formed ahead of the protuberance causes 
the boundary layer to separate from the surface upstream of the fin 
location. This results in a separated flow region composed of horse- 
shoe vortices near the surface and a lambda-type shock pattern in the 
plane of symmetry ahead of the fin. The major effects of the shock- 
wave/turbulent boundary- layer interaction are increased pressure loads 
and heat transfer levels on and in the vicinity of the protuberance. 
The existing experimental data suggest that the separated flow is 
very unsteady due to the oscillation of the separation location with 
amplitude larger than the undisturbed boundary layer thickness. This 
relatively large movement of the separation location is attributed to 
the pulsating and scavenging action of the flow comprising the 
multiple horseshoe vortices, results in comparatively large 
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oscillation of the impinging shock wave. In addition to the pressure 
peaks and heating rates at shock impingement on the fin leading edge, 
there are very large pressures and heating rates on the plate surface 
in the immediate vicinity of the fin root. It is believed that such 
shock wave turbulent boundary layer interactions have caused 
structural damage on hypersonic wing/body junctions [15]. 

Certain aspects of the separated flow ahead of a blunt fin are 
somewhat similar to two-dimensional separated flows ahead of ramps or 
steps. The pressure increase across the fin bow shock wave presents a 
strong adverse pressure gradient to the boundary layer flow on the 
flat plate surface, causing the boundary layer to separate adjacent to 
the plate surface. The effects of pressure increase are propagated 
upstream through the subsonic pocket of the boundary layer and leads 
to the separation of the boundary layer on the flat plate. As 
sketched in Fig. 2.2, the reversed flow forms a horseshoe vortex in 
which the flow quickly spirals downstream, away from the plane of 
symmetry. Immediately downstream of this vortex, there appears to be 
one or more additional vortices. The reversed flow in the vortices is 
constantly replenished by the separated stream flow, and spirals down- 
stream very rapidly. Thus, the vortices bring the high energy air 
stream into proximity with the surface. Extremely high heating rates 
and large pressure fluctuations have been observed in local regions on 
the surface adjacent to the protuberance. The existing oil flow 
visualization techniques suggest that the direction of flow in the 
vortices is predominantly outboard and bears little similarity to two- 
dimensional separated flows in the plane of symmetry. The predomin- 
ance of the transverse flow masks the nature of the flow along the 
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Fig. 2.2 Flow structure at the high speed wing/body junction. 
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boundary between two vortices, but there appears a third, small and 
counter rotating horseshoe vortex between the two large vortices. The 
flow at the root of the fin is highly three-dimensional and complex. 
There are large pressure, as well as high velocity gradients in this 
region. The compressed separated flow ahead of the fin accelerates to 
supersonic speeds in escaping around the fin to lower pressure 
regions. A strong vortex starts at the root of the fin and spreads 
out as it follows the fin root downstream. The fin bow shock causes 
the boundary layer to separate from the surface ahead of the fin, 
resulting in a separated flow region composed of horseshoe vortices 
near the surface, and a lambda-type shock pattern ahead of the fin. 

The shock wave emanating from the separated flow region impinges on 
the fin bow shock, and causes intense heating and high pressure local- 
ly around the fin leading edge. 

Young et al. [16], Price and Stallings [17] and Westkaemper [18] 
observed that the extent of turbulent separation is insensitive to 
Mach number, Reynolds number and boundary layer thickness for a wide 
range of flow conditions, for fin diameters that are larger or compar- 
able to the undisturbed boundary layer thickness. Specifically, for 
Mach numbers from 1.2 to 21, turbulent boundary layers separate 
approximately two diameters ahead of the fin leading edge. Laminar or 
transitional boundary layer separation ahead of a fin is much more 
extensive than turbulent separation and depends on both Mach number 
and Reynolds number. Limited data available in the literature indi- 
cate that the extent of separation ahead of a fin increases with 
increasing Reynolds number based on the separation distance. Laminar 
separation lengths exceeding six diameters have been observed by Young 
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et al . [16] for both Mach 3 and Mach 5 flows. However, separation is 
a strong disturbance to a laminar boundary layer and causes earlier 
transition than for attached, undisturbed boundary layers. In virtu- 
ally all of the tests carried out in the earlier study [16], trans- 
ition occurred either ahead of the separation location or in the free 
shear layer bounding the reverse flow region ahead of the fin (trans- 
itional separation). The extent of transitional separation diminished 

toward the turbulent value U_ ♦ 2D) as the separation location 

sep 

approaches the undisturbed boundary layer transitional location of the 
plate boundary layer. Mach number, Reynolds number, and boundary 
layer thickness effects on the extent of laminar or transitional 
separation have not been clearly differentiated. 

Sidney [19] investigated experimentally the supersonic junction 
flow problem using flow visualization technique. The major objective 
of the experiment was to resolve the fine details of the flow struc- 
ture. The variation of primary separation distance was discussed as a 
function of Mach number, Reynolds number and obstacle dimensions. The 
structure of the separated flow upstream of the obstacle was found to 
change with relatively small changes in Reynolds number. Experiments 
have been conducted by previous investigators to study the effect of 
transition properties on the diameter of the fin [17, 18, 20-27], 
sweep angle [17, 28, 29], height of the protuberance [17, 18, 20], 

Mach number [17, 20, 23, 25], and Reynolds number [17, 20, 25, 28, 

29]. Attempts have been made to correlate the experimental data, but 
confusion over the conditions necessary for the model to be considered 
semi-infinite have caused some difficulty, particularly with respect 
to the vertical scale of the flow field. However most of the 
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experiments revealed clearly that the dominant parameters controlling 
the centerline flow field are the diameter and sweep angle of the fin. 
Dolling and Bogdonoff [30] studied experimentally blunt fin induced 
shock wave/turbulent boundary layer interaction. The first objective 
of the experiment was to investigate the spanwise development of the 
disturbed flow field and to determine its dependance on the model 
geometry and incoming flow conditions. The second objective of the 
experiment was to determine the vertical extent of the interaction on 
the fin. The results show that, on the test surface near the fin and 
on the fin itself, the leading edge diameter plays a dominant role in 
determining the interactions scale and characteristics. 

Hung and Buning [31] studied the blunt fin induced shock wave/- 
turbulent boundary layer interaction problem using the 3D compressible 
Navier-Stokes equations. It was noted that varying the incoming 
boundary layer thickness by an order of magnitude showed that the size 
of the horseshoe vortex and the spatial extent of the interaction are 
dominated by inviscid flow and only weakly dependant on the Reynolds 
number. Hung [32] also noted that by changing the blunt fin to a flat 
faced fin the extent of separation increased drastically, and the main 
horseshoe vortex bifurcated into two vortices rotating in the same 
direction. 

2.2.2 High Speed Controlled Juncture Flows 

Very little work is done in high speed controlled juncture flows. 
Price and Stallings [17] and Stollery [33] studied the effect of fin 
sweep in turbulent boundary layers and noted that the separation of 
the flow ahead of the fin can be reduced drastically by leading edge 
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sweep. Gaitonde and Knight [34] studied numerically the effect of 
suction to control the sharp fin/flat plate junction flow. The compu- 
tation indicated that the particular bleed schedule of five percent 
employed in the study was incapable of controlling or significantly 
modifying the undesirable separation and the large vortical structure. 
It may be possible to signif icantly influence the flow field by minor 
modifications to the configuration such as raising the fin a slight 
distance from the flat plate could provide interesting results. 

2.3 Turbulence Modeling for Wing/Body Juncture Flows 

Complex flows of the type found in wing/body junction are marked- 
ly influenced by three-dimensional ity, pressure gradients, and viscous 
and turbulence effects. An important effect of calculating the flow 
described above is turbulence modeling. The simplest model used for 
the calculation of turbulent shear flows are algebraic eddy viscosity 
models [35]. These models have been employed primarily for calculat- 
ing two-dimensional flows. However, the formulation of eddy viscosity 
models contain a large amount of empiricism that is generally not 
valid in complex three-dimensional flows. Three other types of turbu- 
lence models have been developed in the literature that are applicable 
to junction flow, such as two-equation models [36], algebraic Reynolds 
stress models [37], and the full Reynolds stress transport equations 
[38]. Of the models available, only the algebraic eddy viscosity 
model and the two equation models have been employed in the compu- 
tation of wing/body junction flow [39]. The application of the full 
Reynolds stress transport equations to three-dimensional flows involve 
the solution of six additional partial differential equations for the 
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stresses and greatly increase the computational effort required. The 
algebraic Reynolds stress models have been developed as an attractive 
alternative to the full Reynolds stress equations. However, this 
model has not been refined sufficiently to make them applicable to 
three-dimensional flow fields. In particular, more work needs to be 
done before these models can adequately compute the near wall 
turbulence. 

The literature review indicates that most of the experiments 
carried out in the junction flow lack detailed measurements in the 
separated flow region and do not sufficiently portray the structure of 
the vortical flow. The same level of criticism can be levelled 
against the previous computational studies in that they fail to 
examine the effect of grid refinement and turbulence model in the flow 
field simulation. Although techniques for controlling the junction 
flow have been studied experimentally and computationally in the past, 
quantitative information about their effect on the vortex strength is 
lacking in most of the earlier studies. 

The primary objective of the present study is to obtain further 
understanding of the structure of resulting three-dimensional sepa- 
ration and to compute the junction flow with sufficient resolution. 

The specific objective of the study is to investigate numerically the 
effect of fillets and sweep in controlling the juncture horseshoe 
vortex and to provide quantitative information on the amount of weak- 
ening of the vortex. 



Chapter 3 

PHYSICAL MODEL AND THEORETICAL FORMULATION 
In this chapter, a brief description of the governing partial 
differential equations representing the flow, along with the boundary 
and initial conditions employed in the present study are presented. 
The theoretical formulation of the problem starts with the nondimen- 
sional form of the compressible three-dimensional Navier-Stokes equa- 
tions cast in the generalized body-oriented coordinate system. Brief 
discussions on the turbulence model, mesh system, and thin-layer 
approximation are also provided. 

3.1 Governing Equations in Non-Dimensional Form 
In the absence of body forces and source terms, the governing 
equations for the conservation of mass, momentum, and energy in a 
stationary control volune v, surrounded by a control surface S, can 
be expressed in integral form as [40] 

J_ fl dv + f (£!i + J»> • H ds = 0 (3.i) 

at J\i Js 


where 


Q = [p » pd, pv, pw, pE] 


E = e. + — (u + v + w ) 
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b = [ b , b_, b 1 
-p -m* -e J 

b = 0 

-p 

-l 

b = p I + M Re c 

-m - » - 

^ = - X div £ I - p [(grad jj) + (grad _u) T ] 

/ v - 1 -l 

b = - y M (Re Pr) y grad e. + pu + M Re f u 

— 0 00 J . 00 — 

The nondimensional i zation of the governing equations is carried out by 
normalizing the cartesian velocity components (u,v,w) by sound speed 
a , the density p by p , the specific internal energy e. and the 

00 OO | 

2 2 

total energy E by a^, and p by P^a^. The viscosity coeffi- 
cients X and u are normalized with respect to the molecular 

viscosity p . The constant y, is the ratio of specific heats, Re 
00 

is the Reynolds number based on freestream velocity, and Pr is the 
Prandtl number. The density p, static pressure p and the internal 
energy e. are related through the following equations of state: 

p = (1/y ) pT (3.2) 

ei = T/y (y-1) 


(3.3) 
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a 2 = T 


(3.4) 


where T is normalized by the freestream static temperature. 

Equation (3.1) is a simple expression for the basic conservation laws 
of fluid flow and is valid for a global control volune and also for 
each locally discretized control volume (i, j, k) shown in Fig. 3.1. 
For the sake of notational convenience, the superscript without a 
prime is used for the transformed coordinate (x x ,x 2 ,x 3 ) interchange- 
able with the local coordinates (c, n, c) and the superscript with a 
prime is used for the global coordinates (x,y,z). A similar notation 
convention is used for all other tensor quantities. For example, 
u 1 » (u 1 .u 2 .u 3 ) ■ (u,v,w). 

The partial differential equation corresponding to Eq. (3.1) 
written in tensor convention is expressed as follows: 


3 r , .i' 2 -. x 3 r x .. 1/2 ii n (3.5) 

— Cl (9) ] + [i _u + b)« g l ] = 0 

St 3X 1 


or 


£, t + 1 F,. = 0 


where 1 = 1 (g) , 1 £ = [ (£ JJ + !)• (g ) ] 


1 / 2 . 


In eq. (3.5) (g 1/2 ) is the Jacobian and £ is the contravari ant 

base vector in the x 1 direction. Here, the Einstein sunmation con- 
vention is used for repeated indices. Since Eq. (3.5) is a 




Fig. 3.1 A hexahedron finite volume v(i,j,k) with coordinate systems 
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homogeneous equation in it can be written as 


q, t + ( 1 ^». = 0* i = I. t ^ 5 


Equation (3.6) provides the basis for the implicit procedure. 
Eq. (3.5) explicitly, one has 


continuity 


1/2 , , 1/2 

[(g) p], t + [pu..(g) £ ] , £ = o 


*j 1 

momentum for x direction 


[ (g) 1/2 pu i ] , t + {p u.% • (g) sf ) + P(9) 1/2 g J (-ji) 


. i * 


1/2 


-i 9u n ' / x 1/2 i 

-M. Re \ (g) 

3 X 0 ' 


M. Re" „ [’V'S) 1 ' 2 9(„') + ^ 19)1,2 = ° 

3x n 3x n 


(3.6) 

Writing 

(3.7a) 


(3.7b) 
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energy 

[(g) 1 / 2 pE], £ + { (pE + p)Mg) 1/2 £ £ 


Y M 96. wo „ M \ i / 2 /» 

—V 1(g) g* - — div u_ u • (g) £* 

Re Prax 1 t 1 ') Re 


. M » w r m ' 3 u m ' / 1/2 £ m . 3u n ' 1/2 . 

[u (g) 9(k') + u 7 ( g ) g (n’)Ji’£ 0 ( 3 - 7c ' 

Re 3x k 3x m ' 

£ £ 
where g^.j is the component of the contravari ant base vector £ 

t 1 

along the x direction and g.. , is the contravariant unit base 
vector. 

The above system of equations is assumed to be valid for turbu- 
lent flows if the molecular transport coefficients are replaced by 


u = 


+ 


(3.8a) 


v 

Pr 



(3.8b) 


where u t represents the turbulent eddy viscosity and Pr fc the turbu- 
lent Prandtl number. The molecular dynamic viscosity is evaluated 
using the Sutherland's law of viscosity. 
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3.2 Turbulence Model 

A turbulence model for the three-dimensional configuration of 
Fig. 2.2 is complicated and is still under development. As z « 
the fin surface influences the turbulence; on the otherhand as z 0 
the turbulence is dominated by the flat plate surface. In the junc- 
tion between the fin and flat plate, the effect of turbulence is truly 
three-dimensional. The two-layer turbulent eddy viscosity model of 
Baldwin and Lomax [35] is extended to generalized coordinates to 
determine the turbulent transport coefficients. This model is 
particularly advantageous for certain shock-induced separated flows 
that contain a multiple region in which the length scales are not 
clearly defined. However, the present work does not account for the 
amplification of turbulent intensity after a sudden compression of the 
flow through the bow shock wave. The implementation of more realistic 
turbulence models for the present work awaits further study. 

In each (J-K) plane, the computational domain is like a corner 
formed by two perpendicular walls as sketched in Fig. 3.2. For the 
sake of convenience coordinates (y,z) are used to define each (J-K) 
plane. The two-1 ayer algebraic turbulence eddy viscosity model is 
applied in each region in the following manner. 

Inner Region 

For the inner region, the eddy viscosity model is expressed as 


e . 

inner 


= p (k l 


2 


D) M 


(3.9) 


Location of F 
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where k = 0.4 is Von Karman's constant, z is Buleev's length scale 
[41-43], D is the Van Driest damping factor and u> is the absolute 

value of the vorticity (u> = |v x v |), i.e., 


0) 




(3.10) 


The general definition of the Buleev length scale is 


*(x,y,z) 



(3.11) 


where s defines the distance from the field point (x,y,z) to a point 
on the boundary, and $ is the pheripheral angle as shown in Fig. 

3.3. For a two-dimensional external corner flow, the Buleev length 
scale i is the distance normal to the surface. For an open corner 
formed by two orthogonal planes, the Buleev length scale is given by 


i 


2yz 

y + z + (y 2 + 



(3.12) 


This formula is designed to account for the size of turbulence eddies 
or the turbulence mixing length near the corner under the influence of 
both walls. 




Outer Region 
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In the outer region the outer eddy viscosity is given by 

G outer = pK ^cp ^wake ^kleb (3.13) 

where K is the Clauser constant, C is an additional constant. 

cp 

The outer wake function, F , , is defined by 

wake 

( i F 

F wake = the sma11er of | maX m 2 X * (3.14a) 

KkWe^a*. 

The quantity F is the maximum value of F(y) that occurs in a 
max 

profile and z is the maximum value of z where the maximum value 
max 

of F outer occurs ' ^ relation for F outer is 

F „uter ' * - 0 < 3 - 14b > 

The wake function is a product of a length scale Z and the 
velocity scale F^.^ representative of the energy-containing eddies 
in the outer portion of the Doundary layer. In effect, the distri- 
bution of vorticity is used to determine length scales so that the 
need to find the boundary layer edge is eliminated. The Klebanoff 
intermittency factor is given by: 



F kleb (A) = 


1 + 5.5 


(3.15) 
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where ^ '* s a constant. The quantity U e is the difference 

between maximum and minimum total velocity in the profile. 

U e ■ (Vu* * v2 - ) MX - ( V u2 * V* ♦ ) m1n (3.16) 

The second term in U e is taken to be zero (except in wakes). The 
numerical value of the constants appearing in the foregoing relation 
are given as follows: 

K = 0.0168, C C p = 1.6, C^ e(3 - 0.3 . 

The inner and outer eddy viscosity laws are combined to form the 
composite turbulence eddy viscosity model in the following manner. 
First, profiles of e inner and e outer are determined on each co- 
ordinate line orthogonal to the boundary. The first point nearest to 
the boundary where £ ^ nner exceeds e ou ^ er is denoted the "cross-over 
point." The turbulent viscosity e is equal to e inner f° r all the 
points in the region between the boundary and the crossover point, and 

is equal to e . for all the points above and including the 
outer 

crossover point. 

In the vicinity of the corner, the (J-K) plane is divided into 
four regions I, II, III and IV as shown in Fig. 3.2. In regions I and 
II, D is evaluated at z = 0. The search for F max and the cor- 
responding 4^ proceeds from the wall, either from y = 0 for 
region IV, or from z=0 for region I. The values of F mflX in region 
II and III are treated as constants equal to the value of F^^ at M 
and N, respectively. 
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3.3 Thin Layer Approximation 

Prandtl's classical boundary layer equations are obtained from 
the full Navier-Stokes equations by performing an order of magnitude 
analysis of the terms and neglecting terms of the order of the inverse 
power of the Reynolds number. This procedure leads to (1) neglecting 
the streamwise diffusion compared to normal diffusion; (2) the 
replacement of a normal momentum equation with the assumption of zero 
pressure gradient across the boundary layer thickness. The concept of 
thin-layer approximation to the full Navier-Stokes equation is similar 
to the boundary layer approximation with the only exception of retain- 
ing all three momentum equations. This procedure makes no assumption 
about the pressure variation across the boundary layer thickness. The 
thin-layer approximation is particularly advantageous for computing 
certain high Reynolds number flows and separated turbulent flows 
because it removes the troublesome singularities encountered at the 
separation points and permits straight forward computation of sepa- 
rated flows. As a result of the thin-layer approximation, the govern- 
ing equations become simple and less complicated than those based on 
approximations that attempt to use streamwise gradients parallel to 
the body surface. The considerable simplification in numerical algo- 
rithm and physical concept motivate the straight forward implementa- 
tion of the approximation into computer codes for high Reynolds number 
flows . 

In the present research work this concept is extended in all 
three coordinate directions for a general coordinate system. All the 
viscous terms containing the cross derivatives 9 2 /9x 1 3x J when i * j 
are neglected while retaining the normal derivatives a 2 /ax 1 9x'^ when 
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„ ’ „ n' i 

3 . 3u 3x 

i 9. n 

3x* L 3x* 3X 

It should be noted, that not all the retained terms are larger 
than the neglected terms. Retention of all the second order normal 
derivative terms make the approximation convenient and consistent such 
that all dominant terms are retained for a general coordinate system 
except at the juncture between two walls. Near the juncture between 
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two walls the neglected cross derivatives can be of the same order of 
magnitude as the normal derivatives. But the flow at the juncture 
carries with it a low momentim fluid and hence neglecting the cross 
derivatives will not alter significantly the basic physics of the flow 
field. When not needed (because of the absence of a thin shear 
layer), the viscous terms can be completely neglected in any particu- 
lar direction, and the present algorithm will reduce back to a thin- 
layer approximation for the presence of viscous shear layers in two or 
only one direction. 


3.4 Mesh System 

For simple configurations like swept and sharp fins on a flat 
plate, an algebraic grid generation technique is used to generate the 
mesh. This technique offers an easy control of the grid spacing and 
distribution. A modified elliptic grid generation technique [44] with 
direct control over the grid spacing is used to generate the mesh at a 
filleted fin/flat plate junction. The equation of the super ellipse 

y fi w 

(— ) + ( — ) = 1, is used to generate the grid for circular and 

a b 

continuous fillets in the plane of synmetry. The three-dimensional 
grid for a sharp fin (Fig. 3.4) is obtained by simply stacking 
horizontal grid in the direction normal to the flat plate surface. 

The three-dimensional grid structure for a swept and filleted fin 
(Fig. 3.5 and 3.6) is obtained by simply rotating the grid in the 
plane of symmetry about the center of the fin leading edge. The grid 
is clustered near the fin and plate surfaces using an exponential grid 
stretching formula to provide adequate resolution of the viscous 
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Fig. 3.4 Mesh distribution on sharp fin and flat plate 








PLAT PLATE 


Fig. 3.6 Mesh distribution on filleted fin and flat plate 
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effects. Since the flow is assumed to be steady and symmetrical only 
one-half of the flow is computed. 

3.5 Boundary Conditions 

The physical boundary conditions to be imposed are problem 
dependent. At present, the 3D Navier-Stokes code has been applied to 
blunt fin/flat plate geometries. In this case, a supersonic laminar 
boundary layer develops on the flat plate. The deflection of the fin 
leading edge generates a bow shock wave which interacts with the 
boundary layer on the flat plate. 

The boundary conditions applied may be categorized as follows: 

3.5.1 Inflow Boundary 

The incoming flow is supersonic (except in the subsonic portion 
of the boundary layer) and the inflow boundary is essentially perpen- 
dicular to the oncoming flow. Consequently all flow characteristics 
point from outside towards the computational domain and the appropri- 
ate boundary condition is to specify u, v, w, p and . For the 3-D 
blunt fin, the upstream boundary layer profile is obtained from a two- 
dimensional compressible boundary layer code [45]. 

3.5.2 Outflow Boundary 

The outgoing flow is assumed to be supersonic (except in the 
subsonic portion of the boundary layer), and the outflow boundary is 
essentially perpendicular to the outgoing flow. In this case all the 
flow characteristics point from inside to outside of the computational 
domain, and thus the outflow variables must be determined from the 
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interior flow solutions by zeroth order extrapolation. The fin is 
assumed to be infinite in length and height so that zero-gradient 
boundary conditions are appropriate at the outer boundaries in the 
corresponding directions. The outer boundary of J = J max is 
situated far enough to avoid any influence on the interaction region. 
Theoretically one can prescribe a boundary layer profile on the 
isolated flat plate at each x-location along the outer boundary as 
shown in Fig. 3.7. These profiles are obtained from solving the 
compressible boundary layer equations for flow over a flat plate. 

3.5.3 Solid Boundary 

At a solid wall boundary the following boundary conditions are 
enforced 

7*_n = 0 

(3.18) 

an 

T = T or k — = q 
w 3n 

where v = (u, v, w) is the cartesian velocity vector, p is the 
static pressure, T is the specified wall temperature, q is the 
specified wall heat flux and k is the coefficient of thermal con- 
ductivity. The boundary condition for the static pressure is an 
approximation to the exact boundary condition derivable from the 
normal component of the momentum equations, and has been applied for a 
variety of flows involving shock wave boundary layer interaction. 
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Fig. 3.7 Boundary layer profiles on the flat plate at J = J max 
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3.5.4 Symmetry Boundary 

The fin is at zero angle of attack, and the following boundary 
conditions apply on the plane of symmetry. 

Y.n_ = 0 

= n (3.19) 

dn 

i£- o 

an 

il.o 

3n 

3.6 Initial Conditions and Blunt Body Starting Conditions 
The simplest way to start a computational solution is by assuming 
that all initial conditions are the same as the incoming flow 
solutions. This may result in over compression near the stagnation 
region followed by an over expansion around the shoulder. This may 
cause the numerical scheme to become unstable, thereby resulting in 
negative density. To avoid this, Pulliam and Stegar [46] used a small 
time step during start-up by gradually enforcing the no slip wall 
boundary condition. Rizk, Chaussee, and McRae [47] treated the bow 
shock as a discontinuity and started with an assumed pressure distri- 
bution, in addition to the soft start technique. In the present 
research work, a simple formula is employed to avoid over expansion 
around the shoulder [40]. As depicted in Fig. 3.8 near the nose 
region, the normal velocity component is set to zero on the surface of 
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Fig. 3.8 Region near the nose for setting blunt body starting 
conditions 
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the body at C, and then increased linearly to the incoming flow 
velocity at a predetermined location, say D. The formula for velocity 
at a point between C and D is 


u_ 

P 


' " ( 1 ‘ £/l Cd) 


(u. - _g n )/g nT1 


(3.20) 


where is the fixed length between two points C and 0 and i is 

the variable length from the surface. The density and energy are set 
equal to the incoming flow conditions. This simple device may be 
interpreted as a defect in normal momentum and hence kinetic energy in 
the small region near the blunt nose. This device alleviates the 
nonphysical flow field development and allows large time steps during 
the start up phase. Though different possible boundary conditions 
were discussed the exact boundary conditions used in the present study 
are the same as given in the previous numerical study of supersonic 
junction flow [40]. 


Chapter 4 

METHOD OF SOLUTION 

In this section, the finite volume formulation is developed for 
solving the compressible Navier-Stokes equations at high Reynolds 
number. The governing equations are integrated numerically using the 
MacCormack's [48] explicit-implicit scheme in time split fashion to 
reduce the three-dimensional equations to three sets of locally one- 
dimensional equations. 

4.1 Finite Volume Formulation 

A finite volume cell, indexed by (i,j,k), is defined by eight 

corner points connected with straight lines, as shown in Fig. 3.1. 

Each edge of the finite volume cell is identified by a coordinate 
12 3 

location (x i ,x^ ,x k ) interchangeable with (5-j where 

x{ = ^ = (i-1) Ac, i = 1, IL 

xj = 3 (j-1) An, j = 1 , JL (4.1) 

x k 3 = c k = (k-1) A?, k = 1, KL. 

and the grid spacing along the local coordinate direction is given by 
AX 1 = AC = 1/ (IL-1) 
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Ax 2 = An = 1/ (JL-1) (4.2) 

Ax 3 = A? = 1/ (KL-1) 

The enclosed surface of each finite volume cell consists of a family 
of three coordinate surfaces that delineates the hexahedronal mesh 
cell. The above formula yields a unit cube for the transformed compu- 
tational domain. The choice of a cube of unit dimension in the compu- 
tational domain is arbitrary, and the other investigators have chosen 

the lengths of each side to yield Ax 1 = Ax 2 = ax 3 = 1 so that the 

1/2 

volume of each computational cell is given by dv = g Ax 1 ax 2 ax 3 . 

/ x 1/2 

It should be understood in Eq. (3.5) that g is the volume of the 

1/2 o 

cell dv (i,j,k), the term g g represents the contravari ant 

0 0 

surface normal for the surface x = constant and u . g is the corre- 
sponding contravariant velocity evaluated at the bounding surfaces. 

1/2 o i 

Therefore, in Eqs. (3.7) (g) 9 (i 1 ) 35 the surface area for X = 

] 1 

constant projected on the cartesian coordinate of x = constant. 

, 1/2 1/2 o 

In the finite difference approach, (g ) and g g are 
defined at each grid point location, and are typically evaluated by a 
two-point central difference in all three directions. This leads to 
an inconsistancy in the volirie and surface normal calculations such 
that the geometry conservation law is not satisfied and hence the 
difference scheme cannot recapture the free stream. This inconsis- 
tency does not occur in the present approach. 

Note that any open-surface element for a given boundary enclosing 
the surface has a unique, effective surface vector that is inde- 
pendent of the shape of the surface. This is because, applying the 
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divergence theorem to a constant vector, the integral of the outward- 
oriented surface normal over a closed surface vanishes as given below 

/ * £ ds = / v 7 * d v = 0 for <t> = constant ( 4 . 3 ) 

For example, the surface vector _S 5876 in Fig. 3.1 is independent of 
the choice of which the partitioning surface diagonal is used to 
define cell volumes. Indeed whether all the four vertices lie in a 
plane or not, the effective surface vector is equal to one-half the 
cross-product of its diagonal segments, i.e., 

.§5876 s 0.5 x (r 75 x r 86 ) ( 4 . 4 ) 

Given eight corner points of each finite volume cell, a simple 
way to calculate the volume is to divide each face into two planar 
triangles. The volume is then dependent on which diagonal is used on 
each face since the diagonals of four nonplanar points do not inter- 
sect. In order for neighboring cells to be contiguous, without over- 
laps, neighboring cells must have the same surface partitioning. In 
the present study a simple formula developed by Kordulla and Vinokur 
[ 49 ] is employed to calculate the volume of each cell as 

6 VOL = r 7 i« [ (r 3 1 x £21) + ([21 x [si^ + ([41 * [31) 

+ (r 6 i x r 41 ) + (r 51 x r 61 ) + (r 81 x r 51 )] 

= 1 • [([31 X [2 4 ) + ( ^8 1 X [52) + ([61 X [45) 

= 2 r 71 . (§J432 + §1285 + .§ 1465 ) 


( 4 . 5 ) 
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With the volume and surface evaluated in integral form, Eq. (3.5) 
can be discretized as 


dv. 


i ,J,k 


A q . . , 
i » J , k 


at 


+ 6[(q u + b) . s] 1J#k = 0 (4.6) 


where q. . . is a vector of the dependent variable located at the 

' 9 J 9 K 

center of the cell and (££ + ]>)• .S is the corresponding flux evalu- 
ated at the surface of each finite volune cell. In the present study, 
the discretized finite volume formulation, Eq. (4.6) is used for 
computation. In later discussions, for the purpose of explanation, 

Eq. (3.5) will be used. However, details of terms which have to be 

evaluated are based on Eqs. (3.7a) -(3.7c) with corresponding inter- 

1/2 . .1/2 £ 

pretations of (g) as the volume at the centroid and (g) £ as 

the surface normal evaluated at the boundary of each finite volume 

cell . 

It is interesting to see that the finite volume formulation does 
not require a global coordinate transformation. In fact, the only 
data needed concerning the mesh are the three cartesian coordinates of 
each of the eight vertices of every cell in the mesh system. 

Moreover, the surface and volume of each cell are well defined and 
consistently evaluated, in contrast to a finite difference approach 
where the transformation coefficient and the Jacobian are evaluated 
from surrounding points. 

1/2 

In finite volume formulation, only the cell volume (g) and 

1/2 t 

the surface normal (g )£ are needed for inviscid flow compu- 
tations. To evaluate the viscous terms, the distance between two 
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adjacent mesh cells is needed. In the present approach, instead of 
calculating this distance from the centroid, the vector between two 
adjacent cells (i,j,k) and (i+l,j,k) (see Fig. 4.1) is set in the same 
direction as the surface normal to the surface S^, and is evaluated 
by 





1/2 


1/2 




(4.7) 


When shape variations are large, a more appropriate formula should be 
used for evaluating the metrics as given by the following formula. 




I l-a + ^ -b + - C l (4 ' 8) 


This implies that in a viscous dominated region, for instance near the 
wall, the grid lines would be required to be nearly orthogonal for an 
accurate evaluation of viscous terms. 


4.2 Explicit MacCormack Algorithm 
The MacCormack's explicit algorithm is implemented using a sym- 
metric sequence of time split operators. The concept of splitting, 
generally known as the method of alternating directions, has been most 
popularly used to split the complex operators into a sequence of 
simpler ones. This concept is applied here to reduce the set of 
three-dimensional equations, Eq. (3.5), into three sets of one-dimen- 
sional equations. 

Equation (3.5) can be split into three locally one-dimensional 


(LOD) operators as follows: 
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Lx: = 0 (no semination ini) i = 1, 2, or 3 (4.9) 


Here, each of the split operators contains the spatial derivatives of 
flux in only one direction. Since, this concept is limited to the LOD 
scheme, the summation convention will not be used in this section. 

Let At be the time step, n the time level, and the flux 
across the cell surface x 1 = constant. The LOD predictor-cor- 
rector scheme in a generalized coordinate system can be expressed as 
follows: 


Predictor : 
expl icit 



* - At 


A+Jf_\ 

AX 1 / 


(4.10a) 


update A^ 0 *^ = ^ + a ^ 


Corrector : 
expl icit 



update 


^ n+l _ l /*n ^n+l 

3 . * _ [ a + A a 

2 V 


(4.10b) 


Here, the overbars in the corrector sweeps indicate that the quan- 
tities are determined with updated predictor values. The A + and A 
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indicate forward and backward two-point differences. Although the 
operators are all one-sided by themselves, the resulting algorithm is 
second order accurate in space and time. 

Each of the Lx 1 operators discussed above consists of a 
predictor-corrector set, and requires no block or scalar tridiagonal 
inversion. Because of the noncommutativity of Lx 1 operators, the 
combined numerical scheme is accurate only to first order. To retain 
the second-order accuracy, a symmetric sequence of LOD operators, 
developed by Strang [50] is used. Specifically, the solution vector 
^ is updated in time according to the following relation: 

t n + 2 = Ln (At n ) U (At^) U (2At 5 ) U (At ? ) Ln (At n ) (4.11) 

By this procedure the solution is advanced two time steps from t to 
t + 2at. 

1. Temporal Accuracy 

In order to retain a second order accurate scheme in time for the 
operator sequence of Eq. (4.11), the total time applied to each 
operator must be the same, i.e.. 

At = At = A t = At (4.12) 

£ n 4 

2. Numerical Stability 

Because of the complexity of the full Navier-Stokes equation, it 
is not possible to obtain a closed-form stability expression for the 
MacCormack scheme applied to these equations. However, the linearized 
set of equations, can be subjected to von Neumann type stability 
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analysis. In order to maintain numerical stability, the time for each 
operator must satisfy the following stability condition: 


at c At 


CFL 


At < At 
n n 


CFL 


( 4 . 13 ) 


At < At 
5 ? CFL 

where the maximum allowable time steps At , At and At , 

5 CFL n CFL c CFL 

are given in Appendix A. 

The explicit MacCormack's algorithm of has several advantages. 
First, the method is second-order accurate in time and space. Second, 
it is relatively easy to code the explicit algorithm. Third, the 
explicit algorithm is easily vectorizable. Fourth, the algorithm is 
robust, in the sense that its overall stability is relatively insen- 
sitive to the initial conditions employed. In practical terms, this 

implies that the initial condition need not be close to the final 
solution. In all the computations to date on three-dimensional 
problems, the flow initial conditions are obtained by simply propa- 
gating the prescribed boundary layer at the outer boundary of J = 
J_. v throughout the mesh. 

Computations involving the compressible Navier-Stokes equations 
sometimes "blow up" because of large gradients, for example when the 
flow field contains strong shock waves. In many instances, it is not 
practical to refine the mesh in these regions, particularly if they 
are far removed from the region of interest. For such situations a 
fourth order damping term is added to the scheme to suppress the 
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rurierical oscillations. This smoothing term has a very small 
magnitude except in regions of pressure oscillations where the trunca- 
tion error is already degrading the computed solution, and should be 
nearly zero as steady state is approached. 


4.3 Implicit MacCormack Algorithm 
For computations involving full Navier-Stokes equations, grid 
lines have to be clustered in the direction normal to the wall. In 
such situations, the explicit CFL stability limit becomes very 
restrictive and a large number of iterations are required to reach a 
steady state solution. The stability analysis of the implicit schemes 
shows that there are no restrictions on the time of integration and 
that it is unconditionally stable. 

Bi-diagonal Scheme and Locally One-Dimensional Time Splitting 

The bi-diagonal implicit scheme in generalized coordinates at 
time level n is introduced into explicit predictor and corrector 
steps as follows: 


Predictor : 

explicit 



implicit 


I-At 


* + r t 


AX 


6 ! 


+ 1 = 


A 'PP 

A £ 


(4.14a) 




n+1 


upd ate 


Corrector: 


expl icit 


ttf * 1 - - 4t 


4 '? 


AX 


A .1 ft 1 n+1* \ 

implicit { I + a t -I 1 J 6'q J1 + * = a q" n+ * 


(4.14b) 


AX 


Update : 


^.n + 1 
d 


The matrices | A | 


1 

2 




are matrices with positive eigenvalues. 



(I^l) M, 


where 




/\ 

The M.. are the similarity transformation matrices which diagonalize 
the Eulerian Jacobian ^ in Eq. (3.6) by setting y 3 x = 0. The 
diagonal terms of the matrix ^pl contains three distinct repre- 
sentative eigenvalues defined by 
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. . 2 uyM „ i i Ax 

= max ( | u. i 1 | + 1 7 - » 0 

RePr pAx 1 At 


i = max 
A u+a 


i i -i i / ^ 

| jj . _g_ ± a (g ) | + 


i i .1/2 , 2uyM g 11 AX 1 


r 


- — , o 


RePr pAx At 


where g” = . The eigenvalues of the diagonal matrix ^ D | 

for i = 1 , 2 , 3 can be expressed as follows: 

1111 
for i=l \\ = X 3 = *4 = * u 


V Va* and lx 5 = lx u-a 


2 2 2 
c \ - tm \ S \ 

2 A 4 


for i = 2 Xj = "X„ = “X 


% = 2> u+a and \ - 2x u-a 


3 3 3 3 

for i=3 X^ - a 3 = A u 


3a 4 = 3 x u+a and A 5 A u-a 


Each transformation matrix of M i for i = (1, 2, 3) is split into 
three matrices which can be expressed as follows: 


M i * M ic M it M 
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where 


1 

0 

0 

0 

0 

-up -1 

P- 1 

0 

0 

0 

-vp -1 

0 

P- 1 

0 

0 

-Wp -1 

0 

0 

p- 1 

0 

a0 

-0u 

-0 v 

-0 w 

0 

a = 

0.5(u 2 +v 2 +w 2 ) and 

0 = 

(y-1) 


1 0 0 

0 pa 0 

0 0 1 

0 0 0 

0 -pa 0 


0 e j e ( 2 «) e (3') u 

0 " el ( 2 • ) e l (1 .) 0 0 


0 " el (3 1 ) 0 


• (r) o 

o 1 
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for i = 2 
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for i = 3 



1/2 

Here, e (^*) = g A ^ j / (g**) is the component of the normal 

contravari ant vector in the i' direction. The M matrices transform 
<5q from conservative to nonconservative variables; the matrices 
account for the orientation of the cell faces, x 1 = constant, 
and matrices then transform front nonconservative to character- 

istic forms. It is interesting to see that the matrices M and M ic 
are the same in the Cartesian coordinate system and non-orthogonal 
curvilinear coordinate system; the matrices account for the 

generalized coordinate transformation. The characteristic variables 
for i = 1, 2, 3 are expressed as follows: 
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for i 3 1 



for i = 2 


6p - (1/a 2 )6 p 

6p + pa6u_ . ^/(g 11 ) 

e ‘ (1 • ) SV ' el (2') Su 
e '(l') 6w ‘ e> (3') <SU 
6p - pa6u . gVg 11 ) 


_ 5p - (1/a 2 ) <sp 

e 2 , 6u - e 2 , 6v 

(2 ) (1 ) 

6p + pa6^ . g 2 /(g 22 ) 

e V) { "- e2 (3') * 

6p - pa6^ . £ 2 /( g 22 ) 


for i = 3 


/\ 
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3C 


6p - (1/a 2 ) 6 p 

e3 (3') {U ' e3 (l') SW 

« 3 (3 1 ) 5V ' e3 (2 ') ** 

1/2 

6p + pa6u . g 3 /(g 33 ) 

1/2 

6 p - pa6u . g 3 /(g 33 ) 


Now it is clear that M. transforms the inviscid portion of F.j 
into characteristic form 69^ only when i = j. This implies that, 
in contrast to the unsplit formulation, LOO time splitting makes 
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variables consistent with the characteristic variables by means of the 
similarity transformation employed in the implicit steps [Eqs. (4.14a) 
and (4.14b)]. 

Two points should be noted here. First, it is an insurmountable 
task to develop the similarity transformation for each entire viscous 

Jacobian To keep the right hand side of Eqs. (4.14a) and 

(4.14b) as simple and straightforward as possible, we follow here the 
suggestion given by MacCormack to consider only the Eulerian Jacobian 
of ^ with the addition of a small number of dominant viscous terms to 
the Eulerian eigenvalues. Second, the implicit procedure is skipped 
whenever the explicit stability conditions are satisfied locally, and 
the algorithm requires virtually no more computation time than an 
explicit method. 

Each split operator described above consists of an explicit- 
implicit predictor-corrector set, is second-order accurate in time and 
space, is unconditional ly stable, and requires no olock or scalar 
tridiagonal inversion. Because of the noncommutativity of split 
operators, the simple combined numerical scheme. Lx 1 Lx 2 Lx 3 is only 
first order accurate. To maintain, the second order accuracy of the 
method, a symmetric sequence of the LOO operators, developed by Strang 
[50] are employed. Here, the solution is advanced from time t to t 
+ 2a t as follows 

3 n+2 = Ln (At) U(At) Lc(2At) Le (At) Ln(At) 


(4.15) 



Each locally one-dimensional split operator is applied line by line in 
each sweep direction and, hence, only needs one level plus three 
single arrays of dependent variables. Without this splitting, the 
MacCormack scheme requires, in general, two time levels of variable 
plus a level of right-hand-side residuals for a total of three levels 
of variables stored. 


4.4 Operation of Numerical Code 
4.4.1 General Procedure 

The computation of three-dimensional shock wave boundary layer 
interactions using the numerical code is accomplished in three phases 
as follows: 

1. Mesh Generation and Initial Flow Condition 

Depending on the investigation of the flow past a 
specific configuration, the user generates the grid using a 
separate grid generation code. Once the coordinates of the 
mesh points are known, the volume and surface normal for each 
finite volume cell can be easily computed. The information 
concerning the storage architecture for the physical coordi- 
nates is provided with the numerical code. The initial 
condition, including the upstream (inflow) profiles, and the 

profiles at J= J on the flat plate are provided by the 
max 

user. Depending upon the particular geometry under consider- 
ation the upstream profile may be a uniform flow, or a 
boundary layer flow generated by solving the two-dimensional 
compressible boundary layer equations. The information con- 



cerning the storage for the flow field and surface normals is 
provided in the numerical code. 

2. Computation of the Three-Dimensional Flow Field Using the 
Numerical Code 

The numerical code is used to integrate the governing 
equation in time until a steady state solution is obtained. 
The user of the code must ensure that reasonable criteria for 
the achievement of steady flow have been met. 

3. Graphic Display of Flow Field 

The user must develop his own graphics code to post 
process the output generated by the three-dimensional Navier- 
Stokes code. In this regard, the user may interface his 
graphics code with the graphics library on the host 
computer. 

4.4.2 Requirements on Grid 

The body-oriented generalized coordinate system must provide 
adequate resolution to unravel the basic physical aspects of the flow 
field. A complete set of mesh criteria for judging the acceptability 
of the grid system does not exist at present, except for the process 
of mesh refinement which is generally not feasible for three-dimen- 
sional Navier-Stokes simulations. The determination of the accept- 
ability of a given grid system, therefore, must be based on experience 
and previous numerical computations. 

There are, however, a number of necessary criteria which may be 
identified as follows: 

1. Resolution of Boundary Layers in the Direction Normal to the 
Boundary. 
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The flow field in the direction normal to the boundary 
is generally characterized by the steep flow gradients. This 
is based on previous experience that a minimum of fifteen 
(15) grid points be placed in the boundary layer depending on 
the amount of resolution needed to resolve steep flow gradi- 
ents. For a turbulent flow, the y + of the first mesh point 
should be located less than three wall units from the 
boundary for accurate resolution of the viscous sublayer. 

The distance of the first mesh point from the wall in the 

2 D 

laminar case should be less than or equal to — 

2. Resolution of Boundary Layers in the Cross flow Direction. 

In the numerical simulation of the three-dimensional 
flow field, a typical maximum grid spacing in the cross flow 
direction is provided. This is kept between 0.26 to 0.86, 
where 6 is the boundary layer thickness to accurately 
resolve the spatial features of the flow field. This 
criteria, however, has not been extensively tested and more 
stringent requirements may need to be imposed depending upon 
the flow configuration. 

3. Resolution of the Boundary Layers in the Streamwise Direction 

The streamwise grid spacing depends strongly on the 
character of the flow field. In the vicinity of interaction 
region, the streamwise grid spacing is determined by the size 
of the interaction region. For nominally, two-dimensional 
shock-boundary layer interaction the size is a function of 
pressure rise, Reynolds nimber, and may be a fraction of the 
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boundary layer thickness. For a three-dimensional shock- 
boundary layer interaction, the size of the interaction 
region is larger than a corresponding two-dimensional inter- 
action with the same pressure rise. 

4. Resolution of the Inviscid Flow 

The criteria for grid spacing in the inviscid region 
depends upon the particular flow configuration, and must be 
determined on a case-by-case basis. For example, the exis- 
tence of a bow shock upstream of the blunt leading edge at 
high Mach numbers may impose a stringent local constraint on 
the inviscid grid spacing. 



Chapter 5 

PHYSICAL CONDITIONS OF INVESTIGATION 
For understanding the flow past a wing/body junction, three 
specific simplified geometries were selected. These include: (1) 
blunt swept fin/flat plate junction, (2) sharp fin/flat plate 
junction, and (3) filleted fin/flat plate junction. The special cases 
of blunt and sharp fins are used to obtain specific results for which 
extensive measurements on velocity and pressure distribution, skin 
friction, and surface yaw angles are readily available in the litera- 
ture. In a typical three-dimensional flow past a fin/flat plate 
junction the various parameters controlling the interaction are D 
(diameter of the fin), a (fin sweep angle), 6 (incoming boundary layer 
thickness), (Mach number), and Re^ (unit Reynolds number). The 
parameters controlling the filleted fin-induced interactions are R 
(fillet radius) and 6 (incoming boundary layer thickness). Since 
three dimensional computations are expensive, only limited parametric 
studies were carried out for selected values of sweep angle, fillet 
radius, Mach number, Reynolds number, and Mesh resolution. The range 
of values in each set of parameters is carefully selected to portray 
its dramatic influence in controlling the intersection flow fields. 

The flow field simulation of the selected cases can be easily 
attempted on existing super computers without requiring excessive 
computational resources. 
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5.1 Laminar Flow Past a Swept Fin/Flat Plate Junction-- 
Validation with Experiment and Parametric Study 

5.1.1 Validation with Experiment 

For the specific problem of flow past a blunt fin/flat plate 
junction, comparative results were obtained corresponding to the 
following experimental conditions [16]: 

M . * 2.95, P 0 .* 22261 N/m 2 , T Q „■ 298 K, 

6 = 0.225 cm, D = 1.875 cm, H = 7.5 cm, A = 0 8 

Figure 5.1 shows a simplified blunt fin/flat plate junction and 
various locations along which the computed results are compared with 
the experimental measurements. The experimental investigation 
includes data on pressure distribution along the fin stagnation line 
and selected spanwise locations. To supplement the measurements, oil 
flow visualization technique is also used to reveal the global 
structure of the flow on the fin and flat plate surface. 

5.1.2 Parametric Study 

The parametric studies were conducted to investigate the effect 
of fin sweep and incoming boundary layer thickness. The conditions 
for which specific results were obtained are as follows: 

Effect of Fin Sweep 

M„ * 2.95, p 0 „= 22261 N/m 2 , T Q „■ 298 K 
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6 = 0.225 cm, D = 1.875 cm, H = 7.5 cm, A = 0* to 45* 

Effect of Incoming Boundary Layer Thickness 

M . = 2-95. P 0 2226 N/n. 2 , T Q .= 298 K 

6 = 0.675 cm, D 3 1.875 cm, H = 7.5 cm, A = 0* and 45* only 

5.2 Laminar Flow Past a Blunt Fin - Conditions for Additional 
Results and Parametric Studies 

5.2.1 Conditions for Additional Results 

The following specific conditions were selected to obtain addi- 
tional results and for parametric studies: 

M. = 2.36, p Q 24000 N/m 2 , T Q w 3 294 K 

6 3 0.186 cm, 0 3 0.635 cm, H = 5 cm 

This specific case was selected to validate the code with the experi- 
mental measurements of velocity distribution in the separated region. 

The experimental measurements were obtained by oil flow visualization, 
Schlieren observations, and laser anemometry and the results are 
available in [52]. 

5.2.2 Parametric Studies 

The parametric studies were conducted by using the conditions of 
Sec. 5.2.1 and changing the various parameters. 
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Effect of Grid Refinement 

To investigate the effect of grid refinement, the conditions of 
Sec. 5.2.1 were used and results were obtained for the following three 
grids: 


30 x 50 x 30, 40 x 100 x 30, 40 x 100 x 50 

Effect of Mach Number 

The freestream stagnation temperature was changed corresponding 
to each Mach number while keeping the other conditions constant as 
given in Sec. 5.2.1. The results were obtained for the following 
cases: 


M = 1.2, 

oo 

M = 2.36, 

oo 

M = 3.5, 

00 


T rt = 286 K 
u » 

T_ = 294 K 

(J oo 

T n = 307 K 
0 00 


Effect of Reynolds Number 

The unit Reynolds ntmber of the freestream is altered by changing 
the tunnel stagnation pressure while keeping the other freestream 
conditions constant as given under Sec. 5.2.1. The results were 
obtained for the following cases: 


Re. 


1.25 x 10 6 /m 


PO 


12,000 N/m 2 , 
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p Q m = 24,000 N/m 2 , Re = 2.5 x 10 6 /m 

p Q ^ = 48,000 N/m 2 , Re^ = 5.0 x l(f/m 


5.3 Turbulent Flow Past a Sharp Fin/Plate Plate Junction 

5.3.1 Flow Past a 16* Sharp Fin 

For the case of a 16* sharp fin induced shock wave/boundary layer 
interaction, comparative results were obtained corresponding to the 
following experimental conditions [53]: 

M * 4.0, P n = 1180822 N/m 2 , T n = 336 K 

OO [J 00 U 00 

6 * 0.46 cm, H = 3.7 cm, 0 = 16* 

The experimental study gathered data on wall pressure, surface yaw 
angle, and skin friction distribution at a specified spanwise location 
throughout the interaction region. The experimental study employed 
oil flow visualization to understand the flow structure, Preston tube 
for measuring skin friction, and cobra probes for flow field yaw 
angles. 

5.3.2 Flow Past a 20* Sharp Fin 

The strong interaction generated by a 20* sharp fin on a flat 
plate is investigated corresponding to the following conditions [54]: 
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M = 2.93, P n = 688320 N/m 2 , T n = 260 K 

00 ’ (J 0° U 00 

6 = 1.27 cm, H = 10 cm, 0 = 20* 

Figure 5.2 shows the axial location of experimental surveys where the 
computed pressure distribution is compared with experiment as a 
function of the Z- coordinate normal to the flat plate. The experi- 
ments were carried out using high Reynolds number supersonic tunnel to 
measure the pressure and yaw angle as a function of the coordinate 
normal to the flat plate. 

5.4 Turbulent Flow Past a Swept Fin/Flat Plate Junction 
For the case of swept fin induced shock wave/turbulent boundary 
layer interaction, comparative results were obtained corresponding to 
the following experimental conditions [17]: 

M, = 3.71, p Qoo = 203055 N/m 2 , T Q a> = 339 K 

5=5 cm, D = 5 cm, H = 35 cm, A = 0* to 60° 

5.5 Turbulent Flow Past a Filleted Fin 
5.5.1 Flow Past an Unmodified Junction 

For the case of low Reynolds number flow past a blunt fin/flat 
plate junction, results were obtained corresponding to the following 
experimental conditions [33]: 


p Q 27000 N/m 2 , T Q 312 K 


M. * 2.5, 
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<S = 5 cm, 0 * 2.5 cm, H = 20 cm 

5.5.2 Flow Past a Modified Junction 

For the modified junction, the effects of fillet radius and in- 
coming boundary layer thickness were investigated. 

Effect of Fillet Radius 

Two circular fillets with a radius (R = 5 cm, 9 cm) and a con- 
tinuous fillet were used to modify the junction. The incoming flow 
conditions were kept as in the blunt fin/flat plate junction case. 

Effect of Incoming Boundary Layer Thickness 

Keeping the fillet radius (R a 9 cm) constant, the incoming 
boundary layer thickness was decreased by an order of magnitude from 5 
cm to 0.5 cm to investigate its effect in controlling the juncture 


vortex. 



Chapter 6 

RESULTS AND DISCUSSION 

The entire results obtained for the physical conditions described 
in Chap. 5 are presented and discussed in this chapter. For the 
physical conditions where experimental results were available, the 
numerical results were obtained to validate the nunerical scheme and 
the computer code developed in this study. Extensive results were 
obtained for other physical conditions to study the behavior of three- 
dimensional separation for various high speed juncture flows. The 
effect of different parameters in controlling the intersection flow 
field were also investigated. 

6.1 Laminar Flow Past a Swept Fin 
6.1.1 Computation of 0* Swept Fin Flow Field 

The incoming flow conditions for this case are given in Sec. 
5.1.1. The computations were carried out using a 40 x 40 x 40 mesh. 
Figure 6.1 shows the comparison of pressure distribution along the fin 
leading edge. The leading edge pressure is normalized by the pitot 
tube pressure behind the bow shock wave. The results show that the 
computed pressure distribution agrees well with the experimental 
values. The existence of a peak pressure on the fin leading edge is 
noticed around Z/D * 1.5. The minimum pressure is seen to occur 
around Z/D = 0.25. In the downward direction from this point, the 
pressure is seen to increase slowly towards the plate surface. The 
computed results, therefore, appear to demonstrate the physically 
realistic trend. 
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Fig 6.1 Pressure distribution on the blunt fin leading edge and 
comparison with experiment 
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As pointed out by Hung and Buning [31J, the existence of peak 
pressure on the fin is due to multicompression of the separation 
shock. The decrease of pressure from the peak to the minimum is due 
to extreme expansion created by the horseshoe vortex leading to a 
reversed supersonic zone on the fin. The pressure increase in the 
direction towards the plate surface is due to recompression from 
supersonic to subsonic flow. This creates an adverse pressure 
gradient in the direction toward the plate, resulting in a supersonic 
zone on the fin. 

The computed pressure distribution on the flat plate along the 
line of symmetry is compared with the experimental measurement for an 
unswept fin in Fig. 6.2. The results show good overall agreement 
except near the juncture point. The main features of the flow such as 
upstream influence, pressure plateau, and the pressure rises across 
the fin bow shock and plate separation shock are all well simulated by 
the nimerical code. The existence of low pressure following the sepa- 
ration pressure rise is due to the reversed high speed flow leading to 
a reversed supersonic zone above the plate surface. 

6.1.2 Parametric Effect of Fin Sweep and Boundary Layer Thickness 

The pressure distribution on the fin leading edge is illustrated 
in Fig. 6.3 for various fin sweep angles corresponding to the con- 
ditions given in Sec. 5.1.2. The results are compared with the blunt 
fin case using an identical mesh. It is seen that as the fin sweep 
increases, the location and magnitude of the peak and minimum pres- 
sures decrease drastically. As pointed out by Edney [55], the 
magnitude of the heating rate and peak pressure on the fin leading 





Fig 6.3 Pressure distribution 
fin sweeps 
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edge depend essentially on the geometrical parameters such as the fin 
sweep. As the fin sweep angle increases the supersonic jet becomes 
thinner and may strike the fin leading edge far downstream of the 
impingement point where, because of diffusion its effect on the heat 
transfer and peak pressure may be considerably lower than for the 
unswept case. Edney further observed that as the fin sweep increases 
above a certain critical value, the resulting shear layer may miss the 
fin leading edge entirely, and hence, high pressure and heat transfer 
do not occur. 

Comparison of the pressure distribution on the flat plate is 
shown in Fig. 6.4 for various fin sweep angles. It is noted that as 
the fin sweep angle increases, the separation point on the plate 
surface moves downstream. Also, the peak pressure at the corner 
decreases drastically with increasing sweep angle. Although no 
comparison of pressure distribution has been made with experimental 
measurements for the swept fin case, the computed results agree quali- 
tatively in trend with the Price and Stallings [17] data for a turbu- 
lent boundary layer. 

To further demonstrate the effect of fin sweep on the interaction 
flow field, the code was run for a thick incoming laminar boundary 
layer for conditions given in Sec. 5.1.2. The boundary layer thick- 
ness, in the previous case, was very small. Hence, packing a fine 
grid in this region created some instability in the numerical simu- 
lation. To remedy this situation and also to observe the detailed 
effect of the fin sweep on the interaction flow field, the code was 
run for 0* and 45* fin sweeps. The computations were carried out 
using a 40 x 40 x 40 mesh. 
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Fig. 6.4 Pressure distribution on the flat plate for various 
fin sweeps 
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Figure 6.5 shows the computed pressure distribution on the fin 
leading edge for 0* and 45° swept fins. The results clearly demon- 
strate the dramatic influence of the fin sweep on the pressure distri- 
bution. It is seen that as the fin sweep increases, the peak pressure 
and its location along the fin leading edge decreases dramatically. 

The location and magnitude of minimum pressure on the fin leading edge 
also shows a drastic variation in the 45* sweep case as compared to 
the 0* sweep case. 

Figure 6.6 shows the comparison of pressure distribution on the 
flat plate ahead of the fin leading edge for 0* and 45* swept fins. 

It is seen that as the fin sweep angle increases the point of separa- 
tion of the incoming flow moves closer to the fin leading edge. Also, 
as the sweep increases, no definite first peak in pressure is 
observed. The final pressure jump across the fin leading edge shock 
also shows a drastic reduction in the pressure for the 45* swept fin 
as compared to the 0* swept fin. This is because the flow goes 
through a severe adverse environment, since the fin bow shock is 
aligned normal to the incoming flow for the 0* sweep fin. In the case 
of the 45* swept fin, the incoming flow goes through the oblique bow 
shock formed ahead of the fin and suffers less severe change as 
compared to the 0* swept fin. This is the reason for reduction in the 
level of the peak pressure of the flow ahead of the 45* swept fin. 

Figure 6.7a and 6.7b show the Mach number contours in the plane 
of symmetry for 0* and 45* swept fins. The separation shock, fin bow 
shock, and the two reversed supersonic zones (in the case of 0* sweep) 
have been clearly identified. As the sweep angle increases, the Mach 
number contours show clearly the decrease in the extent of the 



Fig. 6.5 Pressure distribution on the fin leading edge for various 
fin sweeps in a thick laminar boundary layer 
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Fig. 6.6 Pressure distribution on the flat plate for various fin 
sweeps in a thick laminar boundary layer 
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separated region. It is noted that the two reversed supersonic zones 
(one on the fin and one on the plate) are eliminated in the 45° swept 
fin case. A large variation in the Mach number is noted behind the 
fin bow shock for the 0’ sweep fin as compared to the 45° swept fin. 

The normalized static pressure contours in the plane of symmetry 
are shown in Figs. 6.8a and 6.8b for 0* and 45* swept fins. The 
regions of high pressure at the fin leading edge and low pressure at 
the corner are clearly seen in the contour plots. As compared to the 
unswept fin case, the pressure behind the shock for the 45° swept fin 
does not show large differences in the peak pressure. These two argu- 
ments demonstrate clearly the ameliorating effect of the fin sweep 
because the leading edge becomes more and more like a slender ellipse 
in the freestream direction with the increase in the sweep angle. 

This means with a swept blunt fin on a flat plate, there is a reduc- 
tion in the upstream influence of the bow shock, and the effect of 
bluntness associated with the fin is compressed into a small inter- 
action region. In fact, the sweep clearly reduces the magnitudes of 
the interaction and hence the pressure levels throughout the whole 

flow field. The corresponding density contours in the plane of syn- 

% 

metry are shown in Figs. 6.9a and 6.9b. 

6.2 Laminar Flow Past a Blunt Fin 

In this section, numerical results are compared with available 
experimental data. Also, the results of grid refinement and other 
parametric studies are presented and discussed. 









.2 
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6.2.1 Comparison of Computations with Experimental Data 

The base computations for this case have been carried out using a 
30 x 50 x 30 mesh. The tunnel flow conditions are given in Sec. 

5.2.1. The computed pressure distribution along the axis of symmetry 
and comparison with the experiment are shown in Fig. 6.10. It is 
apparent that in the experiment the pressure begins to rise at about 
ten diameters upstream of the fin leading edge. The separation 
occurs, in the experiment, further downstream (X/D s - 6.8) between 
the regions of initial pressure rise and pressure plateau. The pres- 
sure reaches a valley downstream of the plateau region and then 
increases across the detached bow shock. The secondary separation in 
the reversed flow occurs between the pressure valley and the pressure 
plateau. The computed results using the base grid (30 x 50 x 30) are 
in reasonable agreement with the experiment in the region upstream of 
the pressure plateau. As observed in the experiment, the location of 
separation and initial pressure rise are closely predicted by the 
computations. However, the agreement between the computed and experi- 
mental results is not very good downstream of the plateau. Further, 
the computed results fail to predict the length of the plateau region 
accurately. 

The computed pressure distribution along the fin stagnation line 
(4> s 0°), and along various $ = constant and X/D = constant lines on 
the fin surface, are shown in Figs. 6.11a and 6.11b. For = 0° and 
45* , the fin leading edge pressure is referenced to the total pres- 
sure, p t 2 » behind the normal shock. Otherwise, it is referenced to 
the freestream static pressure p^. It is clear from the computed 
results that the peak pressure on the fin surface decreases 
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substantially as the flow expands around the blunt fin. The peak 
pressure on the fin surface at <f> = O’ is about 1.1 and at $ = 45° it 
is approximately 0.6. The magnitude of the minimum pressure also 
drops at $ = 45° but not as much as the peak pressure. The pressure 
drops below the value of the total pressure p t2 behind the normal 

shock. The computed results along the fin side face (Fig. 6.11b) show 

that the pressure drops below the value of freestream static pressure. 
As the flow continues to move downstream the extremely low pressure on 
the fin surface vanishes downstream of the shoulder. 

The computed and the experimentally measured streamwise velocity 
profiles along the line of symmetry are shown in Fig. 6.12. The 
streamwise velocity component is nondimensionalized with reference to 
the freestream velocity U ; U is nominally 560 m/sec. The dashed 
lines along the vertical axis show the displaced origin to avoid over- 
crowding of the figures. The streamwise velocity profiles along X/D = 

-8.0 and -7.2 (Fig. 6.12a) indicate preseparation behavior. The 
velocity profile at X/D * -6.4 indicates that the separation occurs in 
close vicinity of X/D = -6.4. In the experiment, the separation 
occurs around X/D = -6.8. However, it was difficult to make repeat- 
able measurements at that location due to the inherent unsteadiness of 
the flow. Also, due to the small scale of the reversed flow region at 
X/D * -6.4, it was found difficult to detect the reversed flow region 
in the vicinity of the primary separation line. The velocity profile 
at X/D = -5.6 clearly shows the reversed flow region. The velocity 
profiles at X/D = -4.8, -4.0, and -3.2 are shown in Fig. 6.12b. These 
results show that the magnitude and height of the reversed flow region 
increases with increasing X/D. The streamwise velocity profiles 
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Fig. 6.12 Comparison of streamwise velocity profile: (a) along 

X/D = -8.0, -7.2, -6.4 and -5.6 
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Fig. 6.12(continued) (b) along X/D = -4.8, -4.0 and -3.2 
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at X/D = -2.4 and -1.6 are given in Fig. 6.12c. The experimental 
measurements show that the maximum magnitude of the reverse velocities 
at these two locations are approximately equal. At these locations, 
computed results deviate significantly from measurement. The height 
of the reversed flow region decreases with increasing X/D. The 
measured streamwise velocity profiles at X/D = -3.2, -2.4, and -1.6 
display significant kinks between the outer edge of the reversed flow 
region and the edge of the boundary layer. This phenomena is 
attributed to the inherent unsteadiness present in the flow. The 
computed velocity profiles do not display such behavior because the 
flow is assumed to be steady. 

The vertical velocity profiles at X/D = -4.8, -4.0, -3.2, -2.4, 
and -1.6 are shown in Figs. 6.13a and 6.13b. The maximum velocity 
occurs around the outer edge of the boundary layer. The magnitude of 
the velocity then starts decreasing in the z-direction towards the 
freestream. The maximum velocity toward the plate and the height of 
the reversed flow region Fig. 6.12a,b,c increases with increasing 
X/D. 

The postulated mean streamline pattern shown in Fig. 2.2 requires 
positive streamwise velocities close to the flat plate slightly up- 
stream of secondary and tertiary separation lines. The flow structure 
observed in the junction demonstrates that across the boundary layer, 
in the vicinity of secondary and tertiary separation lines, the verti- 
cal velocities remain positive. Moreover, the vortical pattern also 
requires that the vertical velocity profile change its sign more than 
once during a streamwise traverse along various z = constant lines. 
However, the present vertical velocity profiles do not show such 
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Fig, 6.13(continued) (b) along X/D = -2.4 and -1.6 
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behavior. It can be hypothesized that the locations of vortex centers 
are time dependent, but the time averaged surface locations where the 
vortices are shed remain approximately constant. Thus, the streamline 
pattern given in the plane of symmetry can be considered as one of the 
possible streamlines which succeed each other in time. Tobak and 
Peake [56] have shown such a postulated sequence of the cyclic vorti- 
cal structures in low speed laminar flow past a circular cylinder. 

The unsteady nature of high speed laminar flow can be similar to low 
speed laminar flow. However, further study is required to understand 
the unsteady nature of this type of flow field. The limiting stream- 
line patterns on the plate surface shown in Fig. 6.14 depicts the 
separation and reattachment lines on the flat plate. The primary 
separation, corresponding to the most upstream separation location of 
the plate boundary layer, occurs roughly at about X/D = -6.2. The 
reattachment line of the flow is located very close to the fin. 

6.2.2 Effect of Grid and Flow Parameters 

In this subsection, the results of the parametric study of grid 
resolution and variation in Mach and unit Reynolds numbers are pre- 
sented. For various cases, the flow field is investigated for the 
freestream conditions given in Sec. 5.2.2. 

The effects of grid refinement on the wall pressure distribution 
along the line of symmetry are shown in Fig. 6.10. These results 
demonstrate that computations using three different grids yield iden- 
tical results except for the variation in the magnitude of the pres- 
sure valley. The computed pressure distribution along $ s 0* and 45* 
are shown in Fig. 6.11a. Further plots of pressure distribution along 
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4> = 90* and downstream of the fin shoulder (X/D = 2.4 and 4.5) are 
shown in Fig. 6.11b. The pressure distribution along <j> = 0’ is most 
sensitive to the grid refinement. 

The sensitivity of the grid on streamwise velocity profiles at 
various locations upstream of the fin leading edge are shown in Fig. 
6.12. The most upstream separation location of the flat plate 
boundary layer is insensitive to the grid refinement. The computed 
results generally appear to move towards the experimental data with 
refined mesh. However, at X/D = -2.4 and -1.6 (Fig. 6.12c), the 
numerical prediction using fine mesh (40 x 100 x 50) still display 
large discrepancies with the experiment. The vertical velocity pro- 
files and comparison with the experimental data for various streamwise 
locations are shown in Fig. 6.13. The coarse grid (30 x 50 x 30) 
solution at X/D - -4.0 and -3.2 show some discrepancy with the experi- 
ment for z-locations less than 2.5 mm. Computations using the fine 
grid are able to improve the prediction with the experiment at these 
locations. 

The numerically simulated particle paths in the plane of symmetry 
for various grids studied are shown in Fig. 6.15. The solutions show 
the presence of two clockwise rotating primary vortices and a third 
counter clockwise rotating vortex at the corner. These results demon- 
strate that the structure of the flow is relatively insensitive to the 
grid employed. 

The grid refinement study discussed in the previous paragraphs 
established that a 40 x 100 x 30 mesh is sufficient to predict the 
experimental measurements. The effect of Mach number on the pressure 
distribution in the plane of symmetry is shown in Fig. 6.16 using a 40 
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x 100 x 30 mesh. The effect of increasing the Mach number is to 
decrease the length of the most upstream separation distance. Also 
the magnitude of the pressure plateau following separation increases 
with an increase in the Mach number. The pressure distribution along 
the fin stagnation line shown in Fig. 6.17 indicate that as the Mach 
number increases the z-location corresponding to the peak pressure 
decreases. However, the level of peak pressure remains constant for 
M^ = 2.36 and 3.5. The minimum pressure on the fin also shows a 
drastic reduction with increase in the Mach number. Also, the vari- 
ation in pressure from a peak to a minimum increases with increasing 
Mach number. 

The streamwise velocity profiles for various Mach numbers are 
shown in Figs. 6.18a and 6.18b at selected X/D locations. The effect 
of increasing the Mach number is to increase the maximum reverse 
velocity towards the flat plate. The simulated particle paths (in the 
plane of symmetry) at various Mach numbers are shown in Fig. 6.19. 

The results of Figs. 6.19b and 6.19c demonstrate that the flow at 
higher Mach numbers (M^ = 2.36 and 3.5) display no abrupt change in 
the vortical structure. However, at M = 1.2 (Fig. 6.19a), only a 

oo 

single vortex rotating in the clockwise direction can be seen as 
opposed to multi-vortex patterns observed at higher Mach numbers. 
Similar trends in flow patterns were observed by 8aker [3] in low 
speed laminar flow where the number of vortices in the junction 
increases with freestream velocity. 

The effect of unit Reynolds number on the pressure distribution 
along the line of symmetry is shown in Fig. 6.20. The effect of 
increasing the Reynolds number is to move the location of the first 
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Fig. 6.17 Pressure distribution along the fin stagnation line for 
various Mach numbers 
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6.18 Streamwise velocity profiles for various Mach numbers 
(a) along X/D = -4.0 and -3.2 













Fig. 6.19(continued) (c) M,* = 3.5 
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peak in pressure slightly inboard toward the fin leading edge. The 
magnitude of the pressure plateau remains approximately constant. The 
pressure distribution along the fin stagnation line (Fig. 6.21) show 
that as the Reynolds number is increased, the z-location corresponding 
to peak pressure moves slightly upward. The minimum pressure around 
Z/D = 0.5 also shows a substantial variation to the increase in the 
unit Reynolds number. 

The streamwise velocity profile at various X/D locations shown in 
Fig. 6.22 demonstrates that the maximum reverse velocity toward the 
plate increases with an increase in the Reynolds number. Also, the 
height of the reversed flow region above the flat plate decreases with 
an increase in the Reynolds number. The computed particle paths, in 
the plane of symmetry (Fig. 6.23), show the effect of Reynolds number 
on the flow structure. The flow structures at two higher Reynolds 
numbers (Fig. 6.23b and 6.23c) show the presence of three vortices, in 
contrast to the situation at the lowest Reynolds number (Fig. 6.23a) 
which depicts only two vortices. 

6.3 Computation of Three-Dimensional Sharp Fin Flow Field 
6.3.1 Flow Past a 16* Sharp Fin 

The flow past a 16* sharp fin is simulated using the conditions 
given in Sec. 5.3.1. The computations were performed using a 40 x 45 
x 45 mesh. The y + location of the first mesh point ahead of the 
interaction region is set at a distance of less than 5 wall units from 
the fin and flat plate surface. The present computations are compared 
with the previous numerical work of Horstman [56] and the experimental 
measurements of Peake [53] for identical flow conditions. 
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Fig. 6.21 Pressure distribution along the fin stagnation line for 
various unit Reynolds number 
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Fig. 6.23 Particle paths in the plane of symmetry for various unit 
Reynolds number: (a) Re^ = 1.2xl0 6 m _1 




Fig. 6.23(continued) (b) Re^ = 2.5xl0^m"l 
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Fig. 6.23(continued) (c) Re A = 5xl0^m _ l 
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The comparison of present results with previous numerical and 

experimental investigations are given in Figs. 6.24-6.26. Figure 6.24 

shows the wall pressure distribution on the flat plate. The wall 

pressure is nondimensional i zed with respect to the freestream static 

pressure. The abscissa (x-x $hock ) is normalized with respect to the 

incoming boundary layer thickness <5g upstream of the interaction and 

x . , is the axial location of the inviscid shock wave at a specified 

shock 

spanwise station. 

The pressure distribution curve predicts a gradual pressure rise 
due to separation of the incoming boundary layer. The pressure then 
reaches a plateau followed by a final compression as the flow goes 
past the fin leading edge oblique shock wave. All of the features 
observed in the experiment are closely predicted by the present numer- 
ical simulation. It is interesting to see that the present calcu- 
lation agrees very closely with the experiment as compared to the 
previous numerical simulation. A typical inviscid solution is also 
shown in Fig. 6.24. It is clear that the present calculations do not 
give a crisp representation of the shock structure due to the nature 
of the scheme employed in the present study. 

A comparison of the computed surface yaw angles, tan * (v/u), 
with the experimental data is shown in Fig. 6.25. The maximum 
computed surface angle (from the present computation) is about 67.5* 
whereas previous computations have predicted an angle of about 88*. 

The maximum measured angle from the experiment is about 48*. Peake 
[53] observed some nonuniformities in the experiment due to upstream 
nozzle effects. A vortex generator was used upstream of the test 
section to improve the tunnel flow. Because of these conditions, the 
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Fig. 6.24 



Pressure distribution on the flat plate along Y/6 0 = 5.56 
for 16° sharp fin 
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• EXPERIMENT (Ptakc. 19761 



Fig. 6.25 Yaw angle distribution on the flat plate along Y/6 0 = 5.56 
for 16° sharp fin 
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quality of the = 4 experimental data may not be satisfactory. The 
symbols with arrows represent the uncertainty in the experimental 
data. The maximum difference between the computed results and the 
mean experimental data is about 40%. The computed skin friction 
results from the present and previous nunerical simulations are com- 
pared with the experimental data (Fig. 6.26). The experimental 
results are seen to fall between the two sets of computed results. 

6.3.2 Flow Past a 20* Sharp Fin 

For the case of 20° sharp fin induced shock wave/boundary layer 
interaction, comparative results were obtained corresponding to the 
conditions given in Sec. 5.3.2. The computations employed a 30 x 50 x 
50 mesh. The y location of the first mesh point is set at a distance 
of less than 5 wall units from the fin and plate surface. The present 
computational results are compared with the previous computational 
study of Knight et al. [57] and the experimental work of Shapey [54] 
for identical flow conditions. 

The computed pitot pressure profile at station 1 and 2 (Fig. 5.2) 
located upstream of the interaction region are given in Figs. 6.27a 
and 6.27b. The computed and experimental pitot profiles at stations 3 
through 8 (Fig. 5.2) are displayed in Figs. 6.27c-6.27h. The 
abscissa in all figures is the pitot pressure p p normalized by the 
upstream freestream pitot pressure p . The ordinate represents the 

P» 

distance measured from the flat plate surface normalized by the up- 
stream boundary layer thickness 6 q. It should be noted that in 
general 6 Q is not the appropriate vertical scaling factor for this 
complex interaction flow field. As suggested by the experimental data 
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Fig. 6.26 Skin friction distribution on the flat plate along Y/6 0 = 
5.56 for 16° sharp fin 
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Fig. 6.27 



Pitot pressure distribution at different axial locations 
for 20° sharp fin: (a) station 1 
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of McClure and Dolling [58] the appropriate vertical scaling factor is 
1/3 

given by Z Re. /6 n where 6 n i s the experimentally measured boundary 
°0 u 


layer thickness immediately upstream of the shock (without the fin 
effects) at the specified spanwise location. The selection of 6 Q as 
the appropriate scaling factor was motivated by the desire to clearly 
demonstrate the vertical extent of the interaction relative to the 
height of the incoming boundary layer at the location upstream of the 
fin interaction. Similarly, the profiles are shown at selected 
streamwise locations x $ /6 0 , where x $ = x-* shock and * shock is the 
location of the theoretical inviscid shock at a specified spanwise 
location. It should be noted that the experiments by Settles and 
Bogdonoff [59], Dolling and Bodgonoff [60], and Lu and Settles [61] 
indicate that the appropriate scaling is given by X $ = X $ Re fi /6 Q . 

Figure 6.27c shows the comparison of computed pitot pressure 
distribution with the experiment at station 3 which is coincident with 
the line of coalescence. The results from the previous numerical 
study are also included. In general the agreement between the present 
results and the previous experimental and numerical results are quite 
good. Figure 6.27d shows the comparison of the results at station 4, 
located approximately one-third of the distance between the line of 
coalescence and the fin oblique shock wave at the specified spanwise 
location. The computed and experimental profiles display a modest 
"overshoot" outside the boundary layer, associated with the compres- 
sion system ahead of the shock wave. The comparison of computed pres- 
sure profiles with experiment along station 5, located at approxi- 
mately two-thirds of the distance between the line of coalescence and 
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6.27(continued) (c) station 3 



Fig. 6.27(continued) (d) station 4 
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the inviscid shock wave at the specified spanwise position is given in 
Fig. 6.27e. The experimental profiles display a slightly S-shaped 
behavior near the wall, which is less apparent in the computed 
profiles at this location. The overshoot in the pitot pressure is 
more pronounced at this location. In both locations, the computed 
profiles of pitot pressure p p are in reasonable agreement with the 
experimental data, and accurately predict the observed overshoot in 
p p . However, the pressure profiles are smeared near the location of 
Z/6 =1.5 due to the shock capturing nature of the algorithm. 

The calculated and experimental pitot pressure profiles at 
station 6 are shown in Fig. 6.27f. Due to the close proximity of this 
station to the fin leading edge shock wave, uncertainties in pitot 
pressure measurements were noticed in the experimental data outside 
the boundary layer edge. Further improvements in the experimental 
procedure are needed to resolve this issue. But within the boundary 
layer, reasonable agreement is obtained between the experimental and 
the numerical results. The S-shaped nature of the profile is again 
apparent, with reasonable agreement between the experimental and 
numerical results. The computed pitot pressure profiles at stations 7 
and 8 are compared with experimental data in Figs. 6.27g and 6.27h. 

The computed results are in reasonable agreement with the experimental 
data with the discrepancy appearing in the inviscid portion due to the 
shock capturing nature of the algorithm and the proximity of station 7 
to the shock (X $ = 0.6 6 Q ). 

On the basis of the success of validation of the present 
numerical results with the available data, further results were 
obtained to study the details of the instantaneous and limiting 
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6.27(continued) (f) station 6 
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Fig. 6.27(continued) (g) station 7 
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streamline patterns. The instantaneous and limiting streamlines 
clearly depict the complex nature of the interaction flow. Figure 
6.28 shows the computed limiting streamline patterns obtained using 
the Baldwi n-Lomax turbulence model. The line of coalescence and line 
of divergence are clearly observed. These specific nunerically simu- 
lated flow structures are in general agreement with the experimental 
results using a kerosene lampblack visualization technique. 

In particular, the experimental line of coalescence at Y/<S Q = 10 (the 
approximate spanwise extent of the experimental kerosene lampblack 
visualization) is x $ = -5.7 6 Q , whereas the computed lines of coales- 
cence using the Baldwin-Lomax model occur at = -6.45 6 q. The 
computations performed by Knight et al. [57] using the Baldwin-Lomax 
and Jones-Launder models predict the lines of coalescence at X $ = -4.5 

6 and -3.5 6 , respectively. 

0 0 

A series of calculated mean streamlines are displayed in Fig. 
6.29. The vertical scale of the figure has been enlarged for clarity. 
In Fig. 6.29a, the particle traces of 18 streamlines are shown. Six 
streamlines originate from the surface, upstream of the interaction, 
at equal spanwise increments and serve to define the line of coales- 
cence. Six additional streamlines are released immediately above the 
previous six at a height of 0.0048 6 Q . These latter streamlines 
clearly cross the line of separation and continue to roll downstream 
in the clockwise direction (as viewed from the fin leading edge). 
Another set of streamlines is released at a height of 0.52 6 Q from the 
plate surface. These particles clearly display a counter clockwise 
rotational motion of the fluid particles. The streamlines clearly 
display a vortical structure. It is interesting to see that the 



Fig. 6.28 Limiting streamline patterns on the plate surface for 
20° sharp fin 
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particles originating from the higher z-locations are swept beneath 
the particles originating near the surface. Another set of additional 
streamlines originating upstream at a height of Z/6 =1.1 are shown in 
Fig. 6.29b. Unlike the particles originating at lower values of z, 
these particles move approximately parallel to the surface. Also, the 
streamlines in this case exhibit considerable rotation within the 
boundary layer as compared to the case of the 16* sharp fin, due to 
the strong nature of the interaction. Figure 6.30 shows the computed 
limiting streamlines on the fin surface. This again depicts the exis- 
tence of separation on the fin surface. Figure 6.31 shows the 
computed velocity magnitude contours on the flat plate. The simulated 
flow features such as separation shock, fin shock, and reattachment 
shock are clearly evident. 

6.4 Turbulent Flow Past a Swept Blunt Fin 
For the case of the swept fin induced shock wave/boundary layer 
interaction, comparative results were obtained for the specific con- 
ditions of Sec. 5.4. The computations were carried out using a 40 x 
45 x 45 mesh. The pressure distribution along the fin stagnation line 
are compared in Fig. 6.32 for various fin sweep angles. The pressure 
is nondimensional i zed by the total pitot pressure behind the bow 
shock, p t2 » The agreement between the numerical results and the 
experimental data is seen to be good except for the zero sweep blunt 
fin where the level of pressure is underpredicted in the boundary 
layer and overpredicted in the inviscid region. The experimental 
measurements [17] indicated that 
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Fig. 6.30 Limiting streamline patterns on the fin surface for 20 ° 
sharp fin 
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Fig. 6.31 Velocity magnitude contours on the plate surface for 20° 
sharp fin 
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Fig. 6.32 Comparison of pressure distribution along the fin 
stagnation line for various fin sweeps 
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"an irregular pressure distribution does in 
fact occur and the location of peak pressure 
within the boundary layer is a function of sweep 
as well as Mach number. Lack of data in the 
region near the wall prevents identifying a peak 
for any of the other models with this diameter. A 
more complete experimental investigation would 
have to be performed to define, for a certainty, 
the flow model responsible for this irregular 
pressure distribution." 

A close examination of the figure reveals that as the fin sweep from 
the angle increases, the magnitude of the peak pressure decreases 
considerably from 1.2 to 0.2. The variation in pressure from the 
position of peak to minimum is substantially reduced by increasing fin 
sweeps. The increase in pressure in the direction toward the plate 
surface from the position of minimim is not seen for fin sweep angles 
greater than 30*. The computed temperature distribution along the fin 
stagnation line (Fig. 6.33) also shows some improvements due to the 
sweep. 

The comparison of pressure distribution on the flat plate along 
the line of symmetry ( Y/D = 0) and off the line of symmetry (Y/D = 1.5 
to 4.5) are illustrated in Fig. 6.34. Again the agreement is seen to 
be very good. The code is able to simulate all the pertinent flow 
features such as upstream influence, pressure rise due to flow separa- 
tion from the flat plate, pressure plateau and double peaks in pres- 
sure off the line of synmetry. Further examination of the results 
reveal a low pressure region between the pressure plateau and the 
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final compression along the line of symmetry (Fig. 6.34a). The low 
pressure observed is due to the existence of the reversed high speed 
flow region. However, for sweep angles greater than 30° the low pres- 
sure region is not seen. The disappearance of the low pressure trough 
region for high fin sweep angles demonstrates the weakening of the 
main horseshoe vortex with the sweep. Further examination of the 
results illustrates that the sweep considerably decreases the upstream 
influence. The appearance of double peaks in pressure (Figs. 6.34b- 
6.34e) on the flat plate is due to the primary horseshoe vortex moving 
downstream in a helical fashion. The low pressure region between the 
two pressure peaks is caused by a high velocity under the core of the 
vortex [31]. It is seen that the sweep decreases the spanwise extent 
of the interaction considerably. This is due to diminishing of the 
double pressure peaks for higher sweeps. Some of the trends exhibited 
in Fig. 6.34 will become apparent after discussions of the flow struc- 
ture. The computed surface temperature distribution along the line of 
symmetry is shown in Fig. 6.35. The improvement in the surface temp- 
erature is clearly apparent. 

The computed particle paths in the plane of symmetry are shown in 
Fig. 6.36 for blunt and 60* swept fins. The results clearly show the 
flow separation on the flat plate ahead of the fin, leading to the 
formation of the primary horseshoe vortex (Fig. 6.36a). A small 
secondary vortex is also seen near the juncture formed by the fin and 
flat plate (Figs. 6.36a and 6.36b). The secondary vortex rotates in a 
counterclockwise direction as compared with the primary horseshoe 
vortex. The primary vortex elongates to a length of two diameters 
ahead of the blunt fin with its core located about 0.2 diameters above 
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the flat plate. As the fin sweep is increased to 60°, the extent of 
the separation on the flat plate is reduced considerably (Fig. 6.36c 
and Fig. 6.37). This leads to the weakening of the horseshoe vortex 
(Fig. 6.36c). 

The computed skin friction distribution along the line of syn- 
metry is illustrated in Fig. 6.37 for various fin sweeps. It is seen 
that with an increase in the fin sweep, the upstream separation point 
(where value changes sign) moves closer to the fin leading edge. 
Also, for fin sweep angles greater than 0*, no definite reattachment 
point (where the value starts increasing from a negative value) is 
observed. 

The horseshoe vortex transports fresh high momentum fluid contin- 
uously into the separated flow region (Fig. 6.36a). As this high 
momentum fluid accelerates and moves away from the plane of symmetry, 
it leads to the formation of two reversed supersonic zones [31, 62, 
63], one on the fin and the other on the flat plate. The decrease in 
pressure from the peak to a minimum observed in Fig. 6.32 is not only 
due to a lower stagnation pressure in the incoming boundary layer 
profile, but also due to the extreme expansion caused by the horseshoe 
vortex. This expansion leads to the formation of a reversed super- 
sonic zone on the fin. From the computed particle paths in the plane 
of symmetry for the blunt fin (Fig. 6.36a), it is seen that the 
streamlines of low momentum fluid never impinge directly on the fin 
surface. They spiral into the horseshoe vortex and continue to roll 
up downstream in the form of a necklace vortex. The increase in pres- 
sure in the direction toward the plate surface observed in Fig. 6.32 
is due to recompression of the flow from supersonic to subsonic speed. 



158 



Fig. 6.37 Skin friction distribution along the line of symmetry 
for various sweeps 
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It is this adverse pressure gradient that leads to the formation of 
the secondary vortex at the corner (Fig. 6.36b). As the fin sweep 
angle increases, the variation in pressure is reduced substantially. 
This weakens the horseshoe vortex which leads to the disappearance of 
the supersonic zones [62, 63]. 

The computed limiting streamline patterns are shown in Fig. 6.38 
for blunt, 30* and 60* swept fins. The separation and reattachment 
lines are also marked in the figure. The foot trace left by the 
horseshoe vortex as the flow goes around the fin is clearly observed 
in the area bounded by the separation and reattachment line. It is 
apparent that the spanwise influence exerted by the horseshoe vortex 
on the juncture flow is reduced considerably for high sweep angles. 

In particular at a sweep angle of A s 60*, the flow beyond the span- 
wise location, Y/D = 3.0 (Fig. 6.38c) is unaffected by the juncture 
effect. However, in the case of the blunt fin (A = 0*), the extent of 
the disturbed flow persists up to Y/D = 6.0 (Fig. 6.38a). 

6.5 Turbulent Flow Past a Filleted Fin 

For the case of the filleted fin induced shock wave/boundary 
layer interaction, comparative results were obtained for the specific 
conditions given in Sec. 5.5. The computations were carried out using 
a 40 x 45 x 45 mesh. The pressure distribution along the fin stag- 
nation line is compared in Fig. 6.39 for unmodified (blunt fin) and 
modified (filleted) junctions. The fin surface pressure is again 
nondimensional ized by the total pitot pressure behind the bow shock 
(p a 0.12 p^o). The location of the peak pressure corresponds to the 
point of impingement of the stagnation streamline on the fin leading 
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Fig. 6.38 Limiting streamline patterns on the flat plate for various 
fin sweeps: (a) sweep = 0° 


CD| -< 


161 



OIK 



>- IQ 






164 


edge. It is seen that the location of the peak pressure on the fin 
leading edge moves upward for the filleted junction. This is because 
the flow approaching the junction is deflected upward instead of being 
rolled into the horsehose vortex. For the case of the small boundary- 
layer thickness (6/D = 0.2), the results illustrate a considerable 
increase in fin leading edge pressure distribution as compared to the 
unmodified junction. Similar trend in results were observed in blunt 
fin/flat plate studies [31] in which decreasing the incoming boundary 
layer thickness by a factor of five showed a considerable increase in 
the pressure distribution on the fin leading edge. 

The pressure distribution along the line of symmetry are compared 
in Fig. 6.40 for various types of junctions. It is seen that the 
pressure drop, following the pressure plateau, decreases as the fillet 
size increases. This indicates a reduction in speed of the reversed 
flow region and overall acceleration of the juncture flow. 

The computed particle paths in the plane of symmetry are shown in 
Fig. 6.41 for unmodified and modified junctions. The influence of 
filleting on the flow structure is dramatic as compared to unmodified 
junction. The particle paths demonstrate that the horseshoe vortex 
deforms and the flow streamlines within the boundary layer lose much 
of their vortical character with proper filleting. In another study 
[32], in which the leading edge shape was changed from a blunt to a 
square fin, the results showed a considerable increase in the separa- 
tion length. This resulted in a bifurcation of the main horseshoe 
vortex into two vortices rotating in the same direction. For the 
present case, the main horseshoe vortex does not bifurcate even though 
the bluntness of the leading edge is changed by filleting. The 
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Fig. 6.40 Pressure distribution on the flat plate for various 
junctions 
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surface limiting streamlines are shown in Fig. 6.42 for unmodified and 
modified junctions. It is apparent that the separation on the fin 
surface is prevented completely in the case of the modified junction 
as compared to the unmodified junction. The computed instantaneous 
streamlines are shown in Fig. 6.43 for unmodified and modified 
junctions employing circular and continuous fillets. These results 
show the influence of fillets on the helical vortex structure. The 
computed velocity magnitude contours (Fig. 6.44) clearly show that 
with filleting the supersonic zones observed in the case of the blunt 
fin/flat plate junction weaken as the fillet size is increased at the 
juncture. The vorticity contours (Fig. 6.45) demonstrate that the 
magnitude of vorticity in the vortex is reduced by a factor of three 
in the case of the modified junction as compared to the unmodified 
junction. 
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Fig. 6.42 
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Surface limiting streamlines for various junctions: 
(a) unmodified junction (6/D = 2) 
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Fig. 6. 42( continued) (d) circular fillet (6/D = 0.2, R/D = 3 




Fig. 6.43(continued) (b) circular fillet (6/D = 2, R/D = 2) 
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Fig. 6.43(continued) (c) circular fillet (s/D = 2, R/D = 3.5) 
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Fig. 6. 43( continued) (d) circular fillet (6/D = 0.2, R/D 
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ni tude contours in the plane of symmetry 
junctions: (a) unmodified junction (6/D = 2) 
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Fig. 6. 44( continued) (c) circular fillet (6/0 = 2, R/D = 3.5) 
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Fig. 6. 44( continued) (d) circular fillet (6/D = 0.2, R/D = 3.5) 
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Fig. 6.45 Vorticity contours in the plane of symmetry for various 
functions: (a) unmodified junction ( 6/D = 2) 
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Chapter 7 

CONCLUSIONS AND RECOMMENDATIONS 
The supersonic flow past a wing/body junction is simulated 
numerically using a three-dimensional Navier-Stokes code. Results 
have been obtained for simplified geometries such as blunt fin/flat 
plate junction, sharp fin/flat plate junction, and filleted junction. 
Extensive parametric studies of grid resolution, Mach and Reynolds 
number variation have been conducted to understand their effects in 
modifying the interaction flow field. The influence of fillets and 
sweep on the computed flow field have been investigated extensively 
for a wide range of geometrical and flow conditions. In general, the 
computational results show that better agreement with experimental 
results can be obtained by refining the mesh. It has been noted that 
the fin sweep reduces considerably the spanwise extent of the inter- 
action flow field. Some improvements in surface temperature are also 
observed due to sweep. The computed particle paths and limiting 
streamlines for swept junctions demonstrate clearly the dramatic 
influence of fin sweep in weakening the juncture horseshoe vortex. 

The results for sharp fin/flat plate junction show that the structure 
of three-dimensional separation is markedly different from that exhib- 
ited by the blunt fin/flat plate junction. The numerical study of 
filleted junctures clearly show that the flow streamlines lose much of 
their vortical character with proper filleting. It has been demon- 
strated that fillets with a radius of three and one-half times the fin 
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leading edge diameter are required to weaken the usual necklace vortex 
interaction for the Mach number and Reynolds number considered in the 
present study. 

Based upon the present study, the following areas of research may 
be identified for further studies: 

1. The sensitivity of grid refinement on the computed flow 
structure should be investigated using a much finer grid than 
used in the present study to define the minimum number of 
grid points required for an accurate resolution of flow 
features. 

2. The simulated structure of the shock is smeared due to the 
nature of the algorithm applied in the present study. A 
numerical technique, such as a high resolution TVD scheme, 
can be incorporated in the code for providing a crisp repre- 
sentation of the shock. 

3. The flow structure should be examined for physical conditions 
when there are four or more vortices present in the flow 
field. 

4. The effect of changes in the turbulence model, such as a 
nonlinear k- model, can be examined to investigate its 
effect on the flow field simulation. 

5. Since the flow is unsteady, one should investigate the oscil- 
latory behavior of the horseshoe vortex. 

6. It would be desirable to examine the chemistry and real gas 
effects using appropriate mathematical models. 



REFERENCES 


1. Norman, R.S., "On Obstacle Generated Secondary Flows in Laminar 
Boundary Layers and Transition to Turbulence," Ph.D Thesis, 
Illinois Institute of Technology, Chicago, December 1972. 

2. Schwind, R., "The Three-Dimensional Boundary Layer Near a Strut," 
MIT Gas Turbine Laboratory Report, 1962. 

3. Baker, C.J., "The Laminar Horseshoe Vortex," Journal of Fluid 
Mechanics , Vol. 95, Part 2, November 1979, pp. 347-36/. 

4. Kaul, U.K., Kwak, D., and Wagner, C., "A Computational Study of 
Saddle Point Separation and Horseshoe Vortex System," AIAA Paper 
85-0182, January 1985. 

5. Shabaka, I.M.M.A., "Turbulent Flow in a Simulated Wing/Body 
Junction," Ph.D Thesis, Imperial College, London, April 1979. 

6. Shabaka, I.M.M.A. and Bradshaw, P., "Turbulent Flow in an Ideal- 
ized Wing-Body Junction," AIAA Journal, Vol. 19, No. 12, February 

1981, pp. 131-132. 

7. Oguz, E. A., "An Experimental Investigation of the Turbulent Flow 
in the Junction of a Flat Plate and a Body of Constant Thick- 
ness," Ph.D Thesis, Georgia Institute of Technology, 1981. 

8. McMahon, H., Hubbartt, J. and Kubendran, L., "Mean Velocities and 
Reynolds Stresses in a Junction Flow," NASA CR-3605, August 

1982. 

9. Hsing, T.D., and Teng, H.Y., "Experimental Study of the Behavior 
of 3-D Turbulent Boundary Layer in a Simplified Wing/Body 
Junction," AIAA Paper 84-1529, June 1984. 

10. Briley, W.R. and McDonald, H., "Computation of Three-Dimensional 
Horseshoe Vortex Flow Using the Navier-Stokes Equations," Pro- 
ceedings of the Seventh International Conference on Numerical* 
Methods in Fluid Dynamics (W.C. Reynolds and R.W. MacCormack, 
Editors), Springer-ver lag, 1981, pp. 91-98. 

11. Kubendran, L.R. and Harvey, W.D., "Juncture Flow Control Using 
Leading-Edge Fillets," AIAA Paper 85-4097, October 1985. 

12. Mehta, R.D., "Effect of Wing Nose Shape on the Flow in a Wing/ 
Body Junction," Aeronautical Journal, Vol. 88, December 1984, pp. 
456-460. 


190 


13. 


191 

Sung, C.H. and Lin, C.W., "Numerical Investigation on the Effect 
of Fairing on the Vortex Flows Around Airfoil /Flat-Plate 
Juncture," AIAA Paper 88-0615, January 1988. 

14. Kubendran, L., Sung, C.H. and Yang, C.I., "Experiments and Code 
Validation for Juncture Flows," AGARD-CP-437, May 1988, Vol. 2, 

pp. 11.1-11.11. 

15. Korkegi, R.H., "Survey of Viscous Interaction Associated with 
High Mach Number Flight," AIAA Journal, Vol. 9, May 1971, pp. 
771-784. 

16. Young, F.L., Kaufmann, L.G., and Korkegi, R.H., "Experimental 
Investigation of Interaction Between Blunt Fin Shock Waves and 
Adjacent Boundary-Layers at Mach Numbers 3 and 5," Aerospace 
Research Laboratories, ARL 68-0214, December 1968. 

17. Price, A.E. and Stallings, R.L., "Investigation of Turbulent 
Separated Flows in the Vicinity of Fin Type Protuberances at 
Supersonic Mach Numbers," NASA TN 0-3804, February 1967. 

18. Westkaemper, J.C., "Turbulent Boundary Layer Separation Ahead of 
Cylinders," AIAA Journal , Vol. 6, July 1968, pp. 1352-1355. 

19. Sedney, R., "A Survey of the Effects of Small Protuberances on 
Boundary Layer Flows," AIAA Journal, Vol. 11, June 1973, pp. 782- 
792. 

20. Sedney, R. and Kitchens, C.W. Jr., "Separation Ahead of Protuber- 
ances in Supersonic Turbulent Boundary Layers," AIAA Journal , 

Vol. 15, April 1977, pp. 546-552. 

21. Voitenko, D.M., Zubkov, A. I., and Panov, Y.A., "Supersonic Gas 
Flow Past a Cylindrical Obstacle on a Plate," Mekhanika Zhidkosti 
Gaza, Vol. 1, January/February 1966, pp. 121-12IT 


22. Voitenko, D.M., Zubkov, A. I., and Panov, Yu A., "Influence of 
Mach Number on Flow in a Three-Dimensional Separation Region," 
Moscow University Bulletin, Ser. I, Mathematics and Mechanics , 
Vol. 23, No. 2, 1968, pp. 115-118. 

23. Uselton, J.C., "Fin Shock-Boundary Layer Interaction Tests on a 
Flat Plate with Blunted Fins at M = 3 and 5," Arnold Engineering 
Development Center-TR-67-113, June 1967. 

24. Lucas, E.J., "Investigation of Blunt Fin-Induced Flow Separation 
Region on a Flat Plate at Mach Numbers 2.5 to 4.0," Arnold Eng- 
ineering Development Center-TR-70-265, January 1971. 

25. Kaufman, L.G., Korkegi, R.H., and Morton, L.C., "Shock Impinge- 
ment Caused by Boundary Layer Separation Ahead of Blunt Fins," 
Aeronautical Research Lab., ARL 72-0118, August 1972. 



192 

26. Dolling, D.S., Cosad, C.D., and Bogdonoff, S.M., "An Examination 
of Blunt Fin-Induced Flows," AIAA Paper 78-159, January 1978. 

27. Dolling, D.S., Cosad, C.D., and Bogdonoff, S.M., "An Examination 
of Blunt Fin-Induced Shock Wave Turbulent Boundary Layer Inter- 
actions," AIAA Paper 79-0068, January 1979. 

28. Winkelman, A.E., "Experimental Investigations of a Fin Protuber- 
ance Partially Immersed in a Turbulent Boundary Layer at Mach 5," 
(Naval Ordinance Laboratory), NOLTR-72-33, January 1972. 

29. Gillerlain, J.D., "Fine-Cone Interference Flow Field," AIAA Paper 
79-0200, January 1979. 

30. Dolling, D.S. and Bogdonoff, S.M., "Blunt Fin-Induced Shock Wave/ 
Turbulent Boundary-Layer Interaction," AIAA Journal, Vol. 20, No. 
12, December 1982, pp. 1674-1680. 

31. Hung, C.M. and Buning, P.G., "Simulation of Blunt-Fin Induced 
Shock-Wave and Turbulent Boundary-Layer Interaction," Journal of 
Fluid Mechanics , Vol. 154, May 1985, pp. 163-185. 

32. Hung, C.M., "Computation of Separation Ahead of Blunt Fin in 
Supersonic Turbulent Flow," NASA TM-89416, December 1986. 

33. Stollery, J.F., Fomison, N.R., and Hussain, S., "The Effects of 
Sweep and Bluntness on Glancing Interactions at Supersonic 
Speeds," 15th Congress of the International Council of the Aero- 
nautical Sciences, September 1986. 

34. Gaitonde, D. and Knight, D., "Numerical Experiments on the 3-D 
Shock Wave-Boundary Layer Interaction Generated by a Sharp Fin," 
AIAA Paper 88-0309, January 1988. 

35. Baldwin, B.S. and Lomax, H., "Thin Layer Approximation and Alge- 
braic Model for Separated Turbulent Flows," AIAA Paper 78-257, 
January 1978. 

36. Jones, W.P. and Launder, B.E., "The Prediction of Laminarization 
with a Two-Equation Model of Turbulence," International Journal 
of Heat and Mass Transfer , Vol. 15, February 1972, pp. 301-314. 

37. Gorski, J.J., "Computation of the Turbulent Flow in the Corner of 
a Wing/Body Junction," M.S. Thesis, Department of Aerospace Eng- 
ineering, Pennsylvania State University, University Park, 

December 1983. 

38. Hanjalic, K. and Launder, B.E., "A Reynolds Stress Model of Turb- 
ulence and its Applications to Asymmetric Shear Flows," Journal 
of Fluid Mechanics, Vol. 52, April 1972, pp. 609-638. 



193 


39. Gorski, J.J., Govindan, T.R., and Lakshminarayana, B., "Compu- 
tation of Three-Dimensional Turbulent Shear Flows in Corners," 
AIAA Journal , Vol. 23, May 1985, pp. 685-692. 

40. Hung, C.M. and Kordulla, W., "A Time-Split Finite-Volume Algo- 
rithm for Three-Dimensional Flow Field Simulation, AIAA Journal , 
Vol. 22, No. 11, November 1984, pp. 1564-1572. 

41. Hung, C.M. and MacCormack, R.W., "Numerical Solution of Three- 
Dimensional Shock Wave and Turbulent Boundary-Layer Interaction," 
AIAA Journal , Vol. 16, October 1978, pp. 1090-1096. 

42. Buleev, N., "Theoretical Model of the Mechanism of Turbulent 
Exchange in Fluid Flows," Atomic Energy Research Establishment, 
AERE Translation, Hartwell, England, 1963. 

43. Gersner, F., and Po, J., "A Reynolds Stress for Turbulent Corner 

Flows - Part II: Comparisons Between Theory and Experiment, 

Journal of Fluids Engineering, Transaction of ASME, Vol. 98, 
Series I no. 2, June l9/b, pp. 269-277. 

44. Thomas, P.D. and Middlecoff, J.F., "Direct Control of the Grid 
Point Distribution in Meshes Generated by Elliptic Equations," 
AIAA Journal , Vol. 18, No. 6, June 1980, pp. 652-656. 

45. Harris, J.E. and Blanchard D.K., "Computer Program for Solving 
Laminar, Transitional, or Turbulent Compressible Boundary-Layer 
Equations for Two-Dimensional and Axisymmetric Flow," NASA TM- 
83207, February 1980. 

46. Pulliam, T.H. and Steger, J.L., "Implicit Finite-Difference Simu- 
lations of Three-Dimensional Compressible Flow," AIAA Journal 
Vol. 18, February 1980, pp. 159-167. 

47. Rizk, Y.M., Chaussee, D.S., and McRae, D.S., "Numerical Simu- 
lation of Viscous-Inviscid Interactions on Indented Nose Tips," 
AIAA Paper 82-0290, January 1982. 

48. MacCormack, R.W., "A Numerical Method for Solving the Equations 
of Compressible Viscous Flow," AIAA Journal, Vol. 20, September 
1982, pp. 1275-1281. 

49. Kordulla, W. and Vinokur, M., "Efficient Computation of Volume in 
Flow Predictions," AIAA Journal, Vol. 21, June 1983, pp. 917- 
918. 

50. Strang, G., "On the Construction and Comparison of Difference 
Schemes," SIAM Journal of Numerical Analysis, Vol. 5, September 
1968, pp. 555-517; 

Baldwin, B.S., and MacCormack, R., "A Numerical Method for Solv- 
ing the Navier-Stokes Equations With Application to Shock Turbu- 
lent Boundary Layer Interaction," AIAA Paper No. 75-1, January 
1975. 


51. 



194 


52. Ozcan, 0. and Holt, M.," Supersonic Separated Flow Past a Cylin- 

drical Obstacle on a Flat Plate," AIAA Journal, Vol. 22, No. 5, 
May 1984, pp. 611-617. ■ 

53. Peake, D.J., "Three Dimensional Swept Shock/Turbulent Boundary 
Layer Separations with Control by Air Injection," National 
Research Counci 1 -Canada, July 1976. 

54. Shapey, B., "3-D Shock Wave/Turbulent Boundary Layer Interaction 
for a 20° Sharp Fin at Mach 3" MSE Thesis, Mechanical and Aero- 
space Engineering Department, Princeton Unversity, Princeton, New 
Jersey, October 1986. 

55. Edney, B., "Anomalous Heat Transfer and Pressure Distributions on 
Blunt Bodies at Hypersonic Speeds in the Presence of an Impinging 
Shock," Aeronautical Research Institute of Sweden, Stockholm, 
Sweden, Rept. 115, February 1968. 

56. Horstman, C. and Hung, C.M., "Computation of Three Dimensional 
Turbulent Separated Flows at Supersonic Speeds," AIAA Paper No. 
79-0002, January 1979. 

57. Knight, D. Horstman, C., Shapey, B., and Bogdonoff, S.M., "The 
Flowfield Structure of the 3D Shock Wave-Boundary Layer Inter- 
action Generated by 20* Sharp Fin at Mach 3," AIAA Paper 86- 
0343, January 1986. 

58. McClure, W. and Dolling, D., "Flowfield Scaling in Sharp Fin- 
Induced Shock Wave Turbulent Boundary Layer Interaction," AIAA 
Paper 83-1754, July 1983. 

59. Settles, G. and Bogdonoff, S., "Scaling of Two- and Three-Dimen- 
sional Shock/Turbulent Boundary-Layer Interactions at Compression 
Corners," AIAA Journal , Vol. 20, June 1982, pp. 782-789. 

60. Dolling, D. and Bogdonoff, S., "Upstream Influence Scaling of 
Sharp Fin-Induced Shock Wave Turbulent Boundary Layer Inter- 
actions," AIAA Paper 81-0336, January 1983. 

61. Lu, F. and Settles, G., "Conical Similarity of Shock/Boundary 
Layer Interactions Generated by Swept Fins," AIAA Paper 83-1756, 
July 1983. 

62. Lakshmanan, B. and Tiwari, S.N., "Performance of a Three-Dimen- 
sional Navier-Stokes Code on Cyber 205 for High-Speed Juncture 
Flows," ASME Fluids Engineering Conference, Proceedings on Appli- 
cations of Parallel Processing in Fluid Mechanics, FED-Vol. 47, 
June ' 198/, ppT^77T 

63. Lakshmanan, B. and Tiwari, S.N., "Investigation of a Three-Dimen- 
sional Flow Field at a Wing/Body Junction," Deartment of Mechan- 
ical Engineering and Mechanics, College of Engineering and Tech- 
nology, Old Dominion University, Norfolk, Virginia, TR NAS-1- 
17993-03, October 1987. 



APPENDIX A 


MAXIMUM ALLOWABLE TIME STEP FOR EXPLICIT ALGORITHM 


An approximate linearized stability analysis [51] yields the 
following: 


At = CFL .min. At 

t^CFL 1 » * s 

At * CFL .min, At 

n CFL i,J,k n 

(A.l) 


At 
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CFL 4 t { 


where i,j,k denotes the minimum value evaluated over all grid 
volumes which are updated by the explicit scheme, and 
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and 
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u n = ( u n x + v n„ + w n 2 ) / Jn x 2 + n 2 + n 2 2 {A. 6) 

v 1 ".*"'/ * " 'I ' 7 ^X 2 + s * ! z 2 

and 

6 = max |2(vH) - *|, v ( — + — ) (A. 7) 

Pr P rt 

The constant CFL is typically chosen to be 0.9. 
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